Modeling Biology¶

Subcellular Components¶

Daphne is configured with a set of predefined subcellular components: genes, molecules, reactions, and reaction complexes. The user can modify the predefined subcellular components or create new entities.

Genes are defined by assigning a name and specifying the number of copies present in cells and the default transcription activation level (Figure 1a). The activation level can also be dynamically controlled during the experiment, according to the cell state.

Molecules are defined by assigning a name, a value for the diffusion coefficient, and a location, membrane-bound or bulk, (Figure 1b). Bulk molecules can reside in the cytosol or the extracellular medium (ECM), while membrane-bound molecules can only reside in the plasma membrane. By convention, membrane-bound molecules have a | suffix at the end of the molecule name. For example, ‘CXCR4|’ is the name assigned to the membrane-bound CXR4 molecule. In this way, CXCR4 in the cytosol is a different molecule than CXCR4 in the plasma membrane and they can be assigned a different diffusion coefficient than CXCR4 in the cytosol.

Figure 1: Examples of predefined (a) genes and (b) molecules, as viewed in the Subcellular Components Store of the GUI. The user can create new genes or molecules using the ‘Copy’ or ‘New’ buttons or remove items using the ‘Remove’ button

Reactions are defined by specifying a rate constant and selecting the reactant and product molecules that participate in the reaction. Some examples are shown in Figure 2a, including, degradation of CXCR5 molecule , internalization of CXCL13:CXCR5| (the complex formed by CXCL13 binding to CXCR5 receptor) and gene transcription to produce CXCR4 molecule. In this transcription reaction, CXCR4 molecules are created at a rate of , where is the rate constant, is the activity level, andis the copy number.

Users can modify rate constants in predefined reactions and define new reactions (Figure 2b).

Figure 2: .(a) Examples of predefined reactions. The user can edit the rate constants. Together, the first six reactions control homeostasis of CXCR5 receptor. (b) The interface that allows users to create new reactions. To create a new reaction the user selects reactants and products and specifies the rate constant. By definition, any molecule that appears as both a reactant and product is a catalyst.

Currently, Daphne supports twenty reaction types (Figure 3).

Reactions supported by Daphne. In reactions with a mix of upper and lower case letters, a lower case letter represents a bulk molecule and an upper case letter indicates a membrane-bound molecule. In reactions with all lower case letters, the molecules can either be ALL bulk or ALL membrane-bound. The letter ‘g’ denotes a gene

The term Reaction Complex is the name designated for an entity in Daphne that specifies a group of reactions. Typically, the group of reactions in the Reaction Complex represent a set of inter-related reactions that, as a whole, represent a certain process inside the cell or in the extracellular medium. These reaction sets provide a convenient way to add groups of reactions that often occur together to the simulation. Figure 5 shows two examples of predefined Reaction Complexes for CXCR4 and CXCR5 receptor production and recycling that are published with Daphne. In Figure 4, the set of reactions in the selected Reaction Complex, ‘CXCR4 receptor production and recycling’ are listed in the Reactions window. As indicated in Figure 4, users can edit, create, or remove Reaction Complexes by selecting the appropriate button.

Figure 4: Examples of two predefined Reaction Complexes that are published with Daphne

Cell¶

Cells are obviously central to a simulation of the germinal center, and Daphne is published with predefined cell types. The user can modify the predefined cell types or create their own custom cell types. A cell type is uniquely defined by it’s name. In addition to a name, the radius of the cell must be specified, with 5 μm being the default value. In addition to these minimal characteristics, the following optional features can be specified for a cell type: movement, chemistry, and transitions.

Movement¶

Cell locomotion in the simulation is determined by the net force that results from processes such as drag, cell collision, response to chemotactic gradients, cell adhesion, and Brownian motion. The latter is a stochastic process produced by the thermal energy of the surrounding medium.

Brownian (stochastic) motion is parameterized by a parameter that specifies the magnitude the force components

(1)\[ f_i^{(s)} \propto \frac{\sigma}{\sqrt{\delta_t} } z_i\]

where \(z_i\) are independent standard normal deviates and \(\delta_t\) is the time step and the subscript \(i=\) 1,2, or 3 denotes the spatial component of the force vector.

Chemotactic movement is specified by selecting a cytoplasmic driver molecule \(A^*\) and a transduction constant \(\tau\). The chemotactic force experienced by the cell is proportional to the gradient of the driver molecule

(2)\[f^(c) = \tau \nabla A^*\]

where the cytoplasmic driver molecule is denoted by \(A^*\). With appropriate reaction dynamics, inhomogeneity (gradients) in the driver molecule concentration can arise from inhomogeneity in bound chemokine receptor concentrations, which arise due to inhomogeneity in soluble chemokines in the extracellular medium (Figure 5).

Figure 5: Schematic illustration of chemotaxis

For example, the following reaction-based dynamics are used to produce chemotactic motion in many of the predefined cells in Daphne

(3)\[\begin{split}L + R \rightarrow C \\ A + C \rightarrow A^* \\ A^* \rightarrow A\end{split}\]

where \(L\) is a soluble chemokine molecule in the extracellular medium, \(A\) and \(A^*\) are pseudo-molecules in the cytoplasm, and \(R\) and \(C\) chemokine receptor and ligand:receptor complex molecules in the plasma membrane, respectively (e.g., CXCR5 and CXCL13:CXCR5). The first reaction specifies the dynamics for ligand:receptor binding. Inhomogeneity in \(L\) will produce inhomogeneity in \(C\). The second equation specifies the dynamics for activation of pseudo-molecule \(A\) to its activated for \(A^*\). The third reaction ensures that the chemotactic signal will decay in time, as \(A^*\) spontaneously deactivates back to molecule \(A\).

The simple model in (3) is a coarse-grained model that can be replaced by a more detailed intracellular spatial model as the need arises.

Default calibrated parameters are provided for the drag coefficient and cell collision parameters, but the user is also free to override these.

Chemistry¶

The cell genes, molecules, and cell reactions comprise the cell chemistry. Molecules are specified separately for the plasma-membrane and the cytosol (Figure 6).

Figure 6: An example of the GUI for specifying cell molecules and genes.

Cell reactions are, similarly, specified separately for the plasma membrane and the cytosol. Reactions in the membrane can only involve molecules in the membrane. Boundary reactions between plasma membrane and cytosol molecules are, by definition, placed in the cytosol compartment (Figure 7). Similarly, boundary reactions between membrane and ECM molecules are placed in the ECM. Reactions can be added individually or as a group via Reaction Complexes (Figure 7).

Figure 7: Examples of (a) reactions added individually to the Cytosol and (b) reactions add as a group via Reaction Complexes.

Transitions¶

Cell differentiation is specified by defining an Epigenetic Map and Transition Regulators. An example is given in Figure 8, where it can be seen that the (GC B) cell has six possible differentiation states: activated (A), pre-centroblast (pre-CB), centroblast (CB), centrocyte (CC), rescued (R), and apoptotic (Ap). The number and names of the states can be edited by the user, in the case of predefined cells, or defined by the user for custom cells. Figure 8c diagrams the transitions that are possible with the Transition Regulators matrix specified in this example.

Figure 8: An example of a transition scheme for cell differentiation. (a) The epigenetic map, (b) the transition regulators, and (c) a diagram of the possible transitions.

The Epigenetic Map (Figure 9a) allows the user to specify state-dependent activation levels for selected cell genes. In this example, one can see that the gCXCR4 gene is only activated when the cell is in the centroblast state and the gCXCR5 gene is only activated when the cell is in the centrocyte state. The user can add or remove cell genes from the epigenetic map (Figure 9).

Figure 9: Epigenetic map. (a) The Add a gene pulldown menu and (b) user-specified gene activation.

The Transition Regulators table (Figure 9b) allows the user to determine the condition for transition from one state to another. Transitions can be determined by molecular concentrations or probability distributions (Figure 10). Molecule-driven transitions occur stochastically, with increasing probability as the molecular concentration increases. The concentration of the specified driver molecule can be dependent on a single reaction or a complex set of reactions, depending on the biological model that the user wishes to implement. At each time step the probability that the cell transitions from the current state to the destination state is

(4)\[p(t) = \alpha + \beta c(t) \delta t\]

where \(\alpha\) is, background coefficient, and \(\beta\), the linear coefficient, are parameters set by the user, \(\delta t`is the time step, and :math:`c(t)\) is the mean molecular concentration of the selected driver molecule (Figure 10

In the example in Figure 8, the transition from the centroblast to centrocyte state is determined by the molecule W. In this particular example, the concentration of transition driver molecule W starts at a negligible concentration at the beginning of the centroblast state and increases rapidly after a certain number of cell division. Therefore, after a certain number of cell divisions in the centroblast state, there is a large probability that the cell will transition to the centrocyte state.

Figure 10: .(a) A molecule-driven transition. (b0 A probability-driven transition

Also in the example in Figure 8, the transition from the activated to pre-centroblast state is determined by a Weibull probability distribution. That is, the time that a cell spends in the activated state before transitioning to the pre-centroblast state is a number sampled from a Weibull probability distribution. The parameters of the Weibull distribution are set by the user (Figure 10b). Other choices for probability distributions include Gamma, Negative Exponential, Uniform, and Constant (Dirac delta). In the latter, the time to transition is a fixed number specified by the user.

Cell division (the cell cycle) is defined in a similar manner as cell differentiation, with one exception. The last cell cycle state is always defined with a ‘cytokinetic’ state. Cell division occurs when the cell reaches this state, and after cell division the cell automatically returns to the first state. An example of a cell division scheme is shown in Figure 11.

Figure 11: An example of a cell division scheme. Cell division is specified in the same was as cell differentiation, with the exception that the last, cytokinetic, state is fixed. The cell divides upon reaching the cytokinetic state and then returns to the first state (in this example G0) after cell division

Cell death is specified by a single molecule- or distribution-driven event. By default, a dead cell is removed from the simulation in the same time step at which it dies, but the user can specify stochastic removal of dead cells. This can be helpful for calibrating to experimental data.