Because in-person work at a lab was not possible this iGEM cycle, our major advancement was focused on developing a computational model of our system’s circuit design. We wanted our project to delve further into the field, so as a first stride, we created a compartmental, mechanistic ODE-based model using MATLAB©️’s SimBiology. Our model aims to simulate the intracellular reaction dynamics of components of our system, utilizing both mass action and repressor Hill function equation parameters from established biological phenomena, as well as the context of our system--extracellular effects on systemic plasma (Mishra et al.).

Compartmental, mechanistic ODE-based model

Our model is compartmental, in that each physical location of species and reactions is reflected as a space within a membranous organelle or cell compartment. Our model is also mechanistic, in that each step of a molecular process is represented through connections (reaction lines) between input species, binding events, and output species. 

SimBiology is a block-editor based modelling application within the MATLAB language and framework that helps users analyze biological systems characterized by ordinary differential equations (ODEs) (MathWorks 2020). This means that if we can characterize the interplay of species through their rates of change (often well characterized in vitro through kinetics studies in literature) we can make informed predictions on the behavior of a system, as well as iterate critical parameters (ex. rate constants of promoters, leakiness, etc.).

MATLAB Simbiology

We used ODEs to model the production and decay of each component of our system. An example of an ODE in our system is shown below, illustrating the rate of change  of the AND gate Cas6 mRNA transcripts:

$$\frac{d(mRNA_{Cas6})}{dt} = \beta_{CMV} - \delta_{mRNA}[mRNA_{Cas6}] + \frac{\epsilon}{1+\left( \frac{\kappa}{[ERN]}\right )^{n}}$$

where βCMV represents the basal transcription rate of the CMV promoter, δmRNA is the systemic mRNA degradation rate, ε is the ERN’s cut rate, κ is the ERN binding dissociation constant, and n is the Hill coefficient, a parameter crucial to describing the binding behavior of a species.

Intracellular modelling of engineered circuit behavior

If you look into the yellow part of our model figure, this is the circuit architecture that we abstracted in our full circuit diagram. The legend ties this circuit to our original diagram.

Modeling systemic extracellular effects from engineered circuits

To gauge the effect of AND gate output on systemic plasma, we built an extracellular cytokine network map (blue part of our model figure; system output scFv represented as a Y shape, in species blue) based upon principal component analysis (PCA) of cytokine interactions clinically observed in plasma of a population undergoing cytokine storms (Yiu et al., 2012). Rate constants for interactions are coupled concentration coefficients from the study. The interaction-derived  rate constants from literature were used as follows:

From Yiu et al., 2012. Top table shows coefficients of interactions between a pair of cytokines; negative means downregulation, and positive means upregulation. Bottom interaction map shows the three most important upregulators and downregulators (as well as self-attenuation) per cytokine in a cytokine storm model.

This page was written by Sangita Vasikaran, Ethan Levy, Rachel Shen, and Erin Shin