Experiments

Besides modelling the overall system to gain a theoretical understanding as to how the system should behave, we modelled individual experiments so that we could help the biologists and directly compare results by fitting data. We formulated our models using the literature available on the topic as a primary source. This coupled with the simple mass acion law allowed us to formulate models that theoretically should represent the individual experiments quite well. Upon obtaining the results for each experiment, the results would be analysed in the context of the original model's predicted outcome and compared. Due to experimental setbacks this was not possible for all of the experiment models. For some models analytical solutions were obtainable, however for others this was not possible, as a result, we used MATLAB to explore the dynamics of these systems under varying constants. These were the experiments that we modelled:

IRES (Internal Ribosome Entry Site)

We obtained two different IRES sequences whilst constructing our system, one for the NKRF IRES and one for the EMCV IRES. We knew that the EMCV had been shown to work in human cells (Huh-7 cells) [Lohmann 1999] but NKRF had been introduced into HeLa cells successfully [1] and be 92 times more efficient than EMCV in this cell line. We wanted to determine if the same held in Huh-7.5 cells which was our cell line of choice. We modelled the efficiency of the two IRESs and then fitted the data to our model to determine its accuracy. The equations are as follows: \[ \frac{dR_{+}}{dt} = \tau(t) - \mu_{R_{+}}R_{+} \] \[ \frac{dG}{dt} = \beta R_{+} - \mu_{G}G \] This set of equations can be solved analytically to give solutions: \[ R_{+}(t) = \int_{0}^{t} e^{\mu_{R_{+}}(s-t)} \tau(s) ds \] \[ G(t) = \frac{\beta}{\mu_{G} - \mu_{R_{+}}} \int_{0}^{t} ( e^{\mu_{R_{+}}(s-t)} - e^{\mu_{G}(s-t)} ) \tau(s) ds \] where $R_{+}$ is the positive strand of RNA, $G$ is the GFP, $\mu$ is a degradation rate with the subscript denoting the subject, $\beta$ is a translation constant and $\tau$ is the transfection function.

3’ UTR

We found sequences for several different 3’ promoters, and similar to the IRES we wanted to determine which was most efficient. Our goal with the modelling was to determine which would be the best for our system, so that we could then incorporate it into our overall replicon. Due to experimental constraints, the testing of the promoters was only undertaken in E.Coli. The equations governing this experiment are: \[ \frac{dR_{+}}{dt} = c - \mu_{R_{+}}R_{+} \] \[ \frac{dR_{-}}{dt} = \alpha_{-} R_{+} E - \mu_{R_{-}}R_{-} \] \[ \frac{dE}{dt} = 0 \] \[ \frac{dG}{dt} = \beta R_{-} - \mu_{G}G \] where $R_{+}$ is the positive strand of RNA, $c$ is a constant, $\mu$ is a degradation rate with the subscript denoting the subject, $R_{-}$ is the negative strand of RNA, $E$ is RdRp and we assume it's rate of change to be zero indicating that the system has reached equilibrium, $\beta$ is a translation constant and $G$ is the GFP.

Aptazyme

This functions as the kill switch for our system. We wanted to model the effectiveness of this switch, by using values obtained from the literature to obtain theoretical time scales as to how long it would take to stop replication of the system. We then conducted the experiment and fit the data to obtain a comparison as to how close the theoretical model and real values differ. \[ \frac{dR_{+}}{dt} = \tau(t) - \rho X R_{+} - \mu_{R_{+}}R_{+} \] \[ \frac{dX}{dt} = \beta R_{+} - \rho X R_{+} - \mu_{X}X \] where $R_{+}$ is the positive strand of RNA, $\tau$ is the transfection function, $\rho$ is a constant governing the probability of binding, $X$ is the theophylline and$\mu$ is a degradation rate with the subscript denoting the subject.

siRNA

This model was to illustrate the effects of siRNA as a method of gene silencing. We took to the literature and followed previously developed models from these papers, and then explored the dynamics of them using values closer to our systems to gauge a response as to how effective siRNA targeting is at silencing. In the testing this was specific to DPP-IV, however in practice we can take this and apply it to any siRNA targeting sequence. Therefore we made our model general to keep with our aim of reusability so that others in the future could take our model, and with a few tweaks be able to apply it to their systems. \[ \frac{dR_{i}}{dt} = \tau(t) - \rho R_{i} T_{R} - \mu_{R_{i}} R_{i} \] \[ \frac{dR_{ds}}{dt} = \rho R_{i} T_{R} - \lambda_{C} C - ( \lambda_{R_{ds}} + \mu_{R_{ds}} ) R_{ds} \] \[ \frac{dT_{R}}{dt} = - \rho R_{i} T_{R} - \mu_{T_{R}} T_{R} \] \[ \frac{dG}{dt} = \beta T_{R} - \mu_{G} G \] \[ \frac{dC}{dt} = \lambda_{R_{ds}} R_{ds} - \mu_{C} C \] where $R_{i}$ is the "interfering" strand of RNA, $\tau$ is the transfection function, $\rho$ is a constant governing the probability of binding, $T_{R}$ is the target RNA, $\mu$ is a degradation rate with the subscript denoting the subject, $R_{ds}$ is double stranded RNA, $\lambda_{C}$ is a constant, $C$ is the RISC complex, $\beta$ is a translation constant and $G$ is the GFP.

Results

Once the initial models for the experiments had been realised it was clear that some of the dynamics were not being fully captured. Our first ODEs had provided great insight into the way our system operated but the deterministic models did not demonstrate the stability of the system. Simply calculating a trajectory or equilibrium point could easily be very sensitive to initial conditions or changes in parameters. We therefore opted for stochastic DEs. It is a common misconception that the average trajectory of solutions to SDEs will tend towards their ODE counterparts. This will only happen under certain assumptions of stability and so as well as more realistic trajectories, we also got stability analysis.

Next we realised that for any experiment involving RdRp, there would be time when the polymerase was on the mRNA and thus restricting translation. This is particularly relevant to the 3'UTR experiment where the length of the RNA means the polymerase is on for a considerable time. This might be tackled by a new 'combined' RNA and RdRp molecule, however the RdRp is attached to the RNA for a very specific length of time. This lead us to the conclusion that a delay DE would best deal with the dynamics. This is where rate of change is subject not only the current state of the system but also the state of the system in the past. Introducing time delay also brings competition to our model without having to add in probability functions. Our main issue was of combining SDEs and DDEs into one model. Matlab has built in DDE solvers, and there is an well known library that will deal with SDEs. However we could find no m files which would deal with both. So we created it. Our solver runs a modified 4th order Runge-Kutta method which is more accurate than the current tools for SDEs. The dynamics offered by stochastic time delay are much richer than ODEs. These files are available for download .

Our 3'UTR experiment formulated with time delay now looked like; \[ \frac{dR_{+}}{dt} = \tau(t) - \alpha R_{+}(t)E(t) + \alpha R_{+}(t-s)E(t-s) - \mu_{R} R_{+}(t)\] \[ \frac{d[R_{+}E]}{dt} = \alpha R_{+}(t)E(t) - \alpha R_{+}(t-s)E(t-s) \] \[ \frac{dR_{-}}{dt} = \alpha R_{+}(t-s)E(t-s) - \mu_{R_{-}} R(t) \] \[ \frac{dE}{dt} = \beta_{1}R_{+}(t) - \alpha R_{+}(t)E(t) + \alpha R_{+}(t-s)E(t-s) - \mu_{E} E(t) \] \[ \frac{dG}{dt} = \beta_{2}R_{-}(t) - \mu_{G} G(t) \] where the time delays corresponds to the time RdRp takes to copy a strand of RNA.

We then took these equations and adding Weiner noise were able to study how our system behaved in a much more realistic setting. Putting this into our solvers and increasing the number of repetitions to build a more accurate expected value, we found the stochastic model followed the deterministic predictions.