Spatial aspects in biological system simulations

Abstract

Mathematical models of the dynamical properties of biological systems aim to improve our understanding of the studied system with the ultimate goal of being able to predict system responses in the absence of experimentation. Despite the enormous advances that have been made in biological modeling and simulation, the inherently multiscale character of biological systems and the stochasticity of biological processes continue to present significant computational and conceptual challenges. Biological systems often consist of well-organized structural hierarchies, which inevitably lead to multiscale problems. This chapter introduces and discusses the advantages and shortcomings of several simulation methods that are being used by the scientific community to investigate the spatiotemporal properties of model biological systems. We first describe the foundations of the methods and then describe their relevance and possible application areas with illustrative examples from our own research. Possible ways to address the encountered computational difficulties are also discussed.

Figure 17.1

Tabulation of mathematical methods that are often employed in kinetic studies: ODE, ordinary differential equation; PDE, partial differential equation; SPDE, stochastic PDE; SDE, stochastic ODE; MC, Monte Carlo; SSA, Gillespie-type stochastic simulation approach; k_react, rate of biochemical processes; k_diff, rate of diffusion processes.

Figure 17.2

Multicompartment model for receptor signaling. Early endosome (EE) vesicles can form from the cell membrane. EE vesicles can be of two types, smooth and coated pit, with different receptor complex internalization characteristics. Smooth and coated pit EE vesicles also have different trafficking dynamic rates. These vesicles can either recycle back to the plasma membrane or merge into the sorting endosome, which regulates the eventual protein degradation process through the lysosome pathway. The number of EE vesicles changes in time according to the system dynamics. Each compartment is represented with its own receptor signaling network, which is essentially a more complex version of the example given in Eq. (17.1). Further details can be found in Resat et al. (2001, 2003).

Figure 17.3

Schematic diagram showing how spatial details can be incorporated into grid-based multicompartment models during kinetic model development. In this example, cell and its surroundings are modeled as having four compartments: Extracellular (X ), transmembrane (T ), endosomes (E ), and cytoplasm (C ). These subcellular compartments are laid out on a regular rectangular grid where each grid unit is labeled (e.g., E1, X8, etc.) to track location information at the subcompartment level. Note that, if desired, a finer grid may be used for particular regions of the system (e.g., the subdivision of the E1 grid unit).

Figure 17.4

Parsimonious mathematical model to analyze receptor dimerization (Shankaran et al., 2008). The model contains 10 species: Receptor monomers (R₁: EGFR and R₂: HER2) and phosphorylated receptor dimers (R₁₁*, EGFR homodimers; R₁₂*, EGFR-HER2 heterodimer; and R₂₂*, HER2 homodimer) at the cell surface (subscript s) and in the internal (i) compartments. Monomers at the cell surface interact to yield phosphorylated dimers with a forward rate k_fs and a dimer-dependent reverse rate k_rs. Internalized dimers dissociate with a dimer-dependent rate k_ri. Monomers are internalized with rate k_t. Dimers are internalized with a dimer-dependent rate k_e. In the absence of ligand, EGFR and HER2 monomers are assumed to have a surface-to-internal receptor ratios of α₁ and α₂, respectively, and the respective monomer degradation rates k_d1 and k_d2 can be expressed as a product of these α values and k_t. Internalized dimers are degraded with a dimer-dependent rate k_d*. V_R1 and V_R2 are the zero-order synthesis rates for EGFR and HER2.

Figure 17.5

Description of a grid framework in biofilm simulations (adapted from Noguera et al., 1999b). The physical domain is divided into rectangular grid units. First layer forms the surface that microbial species can attach to grow into biofilms. Grid units are labeled as being occupied by the growing biofilm (gray), colony boundary units (yellow). Note that only the units in the first layer are illustrated in the figure for clarity, and the rest are not shown. The diffusing elements can occupy any of the grid units and their dynamics can be modeled using the reaction-diffusion equation, Eq. (17.4). Time evolution of biomass in the grid units is tracked. When a grid unit is filled with biomass, the growing colony occupies the neighboring sites and the interface boundary layer moves accordingly.

Figure 17.6

Construction of multigrid framework in multiresolution coarse-grained approaches. In this particular cell receptor signaling study, three levels of resolution was used: Top lattice (#1) contains seven layers at a z-lattice spacing of 10 nm. First and second layers respectively represent the plasma membrane and part of cytosol that forms the membrane boundary. Next two lattices (#2 and #3) are coarse grained at a z-spacing of 50 nm, and the last lattice (#4) is coarse-grained to 100 nm. Layers 2–7 of lattice 1 and the bottom coarse grained lattices 2–4 represent the cytoplasm. Layers in the lattices use a 50 × 50 grid mesh with 10 nm spacing in the x–y directions, and individual receptors are represented by their occupancy on a single grid unit. When bound to a membrane receptor, cytoplasmic species can diffuse as part of a receptor complex at the first layer of the top lattice. Upon dissociation from the receptor complex, adaptor proteins are placed at the second layer. Cytoplasmic species can diffuse in between lattices as they move around in the cytoplasm.

Figure 17.7

Model for cellulose utilization by a microbial community. Although the actual model is constructed and simulated in 3D, results of a 2D simulation are reported here to keep illustration simple. (A) Metabolic network of the microbial cells. Cells uptake soluble carbon substrate S, and use it for maintenance, protein synthesis, and biomass growth. Cells express (i) regulatory ribosomal proteins EnzP, which controls the growth rates; (ii) transport enzymes EnzT, which facilitates the substrate uptake; and (iii) hydrolases EnzX, which convert polymeric carbon cellulose to simpler carbon forms such as cellobiose, which is consumed by the cells for maintenance and growth. Synthesis rates of these proteins depend on the growth state of the cells. (B) Spatiotemporal profile of cellobiose and biomass in the system. Simulations use a 31 × 31 grid with a 5 micron grid size; the dark blue areas are soil-covered grid points and the light blue units are the pores for the cells to reside. Panels report and compare the soluble substrate and the biomass in the system at the indicated times. Although the simulations keep track of individual cells, total biomass in the grid units is reported using a color scale.