Difference between revisions of "Team:Fudan-TSI/Model"

(Undo revision 494480 by Kiyochou (talk))
Line 1: Line 1:
 
+
{{Fudan-TSI}}
<!DOCTYPE html>
+
<html>
<html lang="en" dir="ltr" class="client-nojs">
+
<head>
+
<meta charset="UTF-8" />
+
<title>Team:Fudan-TSI/Model - 2019.igem.org</title>
+
<meta name="generator" content="MediaWiki 1.24.1" />
+
<link rel="alternate" type="application/x-wiki" title="Edit" href="/wiki/index.php?title=Team:Fudan-TSI/Model&amp;action=edit" />
+
<link rel="edit" title="Edit" href="/wiki/index.php?title=Team:Fudan-TSI/Model&amp;action=edit" />
+
<link rel="shortcut icon" href="/favicon.ico" />
+
<link rel="search" type="application/opensearchdescription+xml" href="/wiki/opensearch_desc.php" title="2019.igem.org (en)" />
+
<link rel="alternate" hreflang="x-default" href="/Team:Fudan-TSI/Model" />
+
<link rel="copyright" href="https://creativecommons.org/licenses/by/3.0/" />
+
<link rel="alternate" type="application/atom+xml" title="2019.igem.org Atom feed" href="/wiki/index.php?title=Special:RecentChanges&amp;feed=atom" />
+
<link rel="stylesheet" href="https://2019.igem.org/wiki/load.php?debug=false&amp;lang=en&amp;modules=mediawiki.legacy.commonPrint%2Cshared%7Cmediawiki.skinning.content.externallinks%7Cmediawiki.skinning.interface%7Cmediawiki.ui.button%7Cskins.igem.styles&amp;only=styles&amp;skin=igem&amp;*" />
+
<!--[if IE 6]><link rel="stylesheet" href="/wiki/skins/Igem/IE60Fixes.css?303" media="screen" /><![endif]-->
+
<!--[if IE 7]><link rel="stylesheet" href="/wiki/skins/Igem/IE70Fixes.css?303" media="screen" /><![endif]--><meta name="ResourceLoaderDynamicStyles" content="" />
+
<style>a:lang(ar),a:lang(kk-arab),a:lang(mzn),a:lang(ps),a:lang(ur){text-decoration:none}
+
/* cache key: 2019_igem_org:resourceloader:filter:minify-css:7:35f114711b15fda1f15f5df02b43e4bc */</style>
+
<script src="https://2019.igem.org/wiki/load.php?debug=false&amp;lang=en&amp;modules=startup&amp;only=scripts&amp;skin=igem&amp;*"></script>
+
<script>if(window.mw){
+
mw.config.set({"wgCanonicalNamespace":"","wgCanonicalSpecialPageName":false,"wgNamespaceNumber":0,"wgPageName":"Team:Fudan-TSI/Model","wgTitle":"Team:Fudan-TSI/Model","wgCurRevisionId":489243,"wgRevisionId":489243,"wgArticleId":6656,"wgIsArticle":true,"wgIsRedirect":false,"wgAction":"view","wgUserName":"Kiyochou","wgUserGroups":["*","user","autoconfirmed"],"wgCategories":[],"wgBreakFrames":false,"wgPageContentLanguage":"en","wgPageContentModel":"wikitext","wgSeparatorTransformTable":["",""],"wgDigitTransformTable":["",""],"wgDefaultDateFormat":"dmy","wgMonthNames":["","January","February","March","April","May","June","July","August","September","October","November","December"],"wgMonthNamesShort":["","Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct","Nov","Dec"],"wgRelevantPageName":"Team:Fudan-TSI/Model","wgUserId":3675,"wgUserEditCount":230,"wgUserRegistration":1562576588000,"wgUserNewMsgRevisionId":null,"wgIsProbablyEditable":true,"wgRestrictionEdit":[],"wgRestrictionMove":[]});
+
}</script><script>if(window.mw){
+
mw.loader.implement("user.options",function($,jQuery){mw.user.options.set({"ccmeonemails":0,"cols":80,"date":"default","diffonly":0,"disablemail":0,"editfont":"default","editondblclick":0,"editsectiononrightclick":0,"enotifminoredits":0,"enotifrevealaddr":0,"enotifusertalkpages":1,"enotifwatchlistpages":1,"extendwatchlist":0,"fancysig":0,"forceeditsummary":0,"gender":"unknown","hideminor":0,"hidepatrolled":0,"imagesize":2,"math":1,"minordefault":0,"newpageshidepatrolled":0,"nickname":"","norollbackdiff":0,"numberheadings":0,"previewonfirst":0,"previewontop":1,"rcdays":7,"rclimit":50,"rows":25,"showhiddencats":0,"shownumberswatching":1,"showtoolbar":1,"skin":"igem","stubthreshold":0,"thumbsize":5,"underline":2,"uselivepreview":0,"usenewrc":1,"watchcreations":1,"watchdefault":1,"watchdeletion":0,"watchlistdays":3,"watchlisthideanons":0,"watchlisthidebots":0,"watchlisthideliu":0,"watchlisthideminor":0,"watchlisthideown":0,"watchlisthidepatrolled":0,"watchmoves":0,"watchrollback":0,
+
"wllimit":250,"useeditwarning":1,"prefershttps":1,"language":"en","variant-gan":"gan","variant-iu":"iu","variant-kk":"kk","variant-ku":"ku","variant-shi":"shi","variant-sr":"sr","variant-tg":"tg","variant-uz":"uz","variant-zh":"zh","searchNs0":true,"searchNs1":false,"searchNs2":false,"searchNs3":false,"searchNs4":false,"searchNs5":false,"searchNs6":false,"searchNs7":false,"searchNs8":false,"searchNs9":false,"searchNs10":false,"searchNs11":false,"searchNs12":false,"searchNs13":false,"searchNs14":false,"searchNs15":false});},{},{});mw.loader.implement("user.tokens",function($,jQuery){mw.user.tokens.set({"editToken":"1d14107e262517d1c37967e6f35eea61+\\","patrolToken":"e787c08cce9ce3fd7234a08755b199b7+\\","watchToken":"e298b55f37d8892d34b476ad3301e9e0+\\"});},{},{});
+
/* cache key: 2019_igem_org:resourceloader:filter:minify-js:7:17aec3a877c68567a6f74963d5548922 */
+
}</script>
+
<script>if(window.mw){
+
mw.loader.load(["mediawiki.page.startup","mediawiki.legacy.wikibits","mediawiki.legacy.ajax"]);
+
}</script>
+
</head>
+
<body class="mediawiki ltr sitedir-ltr ns-0 ns-subject page-Team_Fudan-TSI_Model skin-igem action-view">
+
 
+
        <script type='text/javascript'        src ='/common/tablesorter/jquery.tablesorter.min.js'></script>
+
        <link rel='stylesheet' type='text/css' href='/common/tablesorter/themes/groupparts/style.css' />
+
        <link rel='stylesheet' type='text/css' href='/common/table_styles.css' />
+
 
+
        <script type='text/javascript'        src ='/wiki/skins/Igem/resources/2019_skin.js'></script>
+
 
+
 
+
    <div id='globalWrapper'>
+
        <div id='top_menu_under' class='noprint'></div>
+
        <div id='top_menu_14' class='noprint'>Loading menubar.....</div> <!-- Will be replaced with the jQuery.load -->
+
<script>jQuery('#top_menu_14').load('https://2019.igem.org/cgi/top_menu_14/menubar_reply.cgi',
+
    {  t:"Team%3AFudan-TSI%2FModel",
+
a:"View+%2FTeam%3AFudan-TSI%2FModel++Edit+%2Fwiki%2Findex.php%3Ftitle%3DTeam%3AFudan-TSI%2FModel%26action%3Dedit++History+%2Fwiki%2Findex.php%3Ftitle%3DTeam%3AFudan-TSI%2FModel%26action%3Dhistory++Move+%2FSpecial%3AMovePage%2FTeam%3AFudan-TSI%2FModel++Unwatch+%2Fwiki%2Findex.php%3Ftitle%3DTeam%3AFudan-TSI%2FModel%26action%3Dunwatch%26token%3D997d2d776ae7b3100b2da852a6520e2c%252B%255C++Page+%2FTeam%3AFudan-TSI%2FModel++Discussion+%2Fwiki%2Findex.php%3Ftitle%3DTalk%3ATeam%3AFudan-TSI%2FModel%26action%3Dedit%26redlink%3D1++" });
+
</script>
+
 
+
        <!-- Content div contains HQ_page for HQ styles, Logo and title div, and USER CONTENT -->
+
<div id="content" class="mw-body" role="main">
+
    <a id="top"></a>
+
 
+
            <div id="top_title">
+
                <div class="logo_2019">
+
                    <a href="https://2019.igem.org">
+
                    <img src="https://static.igem.org/mediawiki/2019/8/8b/HQ_page_logo.jpg" width="100px">
+
                    </a>
+
                </div>
+
 
+
        <h1 id="firstHeading" class="firstHeading">
+
            <span dir="auto">Team:Fudan-TSI/Model</span>
+
        </h1>
+
            </div>
+
 
+
            <div id="HQ_page">
+
                <div id="bodyContent">
+
            <div id="mw-content-text" lang="en" dir="ltr" class="mw-content-ltr"><p>
+
  
  
Line 615: Line 550:
 
</head>
 
</head>
  
<link rel="stylesheet" href="https://2019.igem.org/wiki/index.php?title=Template:Fudan-TSI/Fudan-font-awesome.css&amp;action=raw&amp;ctype=text/css" />
+
<link rel="stylesheet" href="https://2019.igem.org/wiki/index.php?title=Template:Fudan-TSI/Fudan-font-awesome.css&action=raw&ctype=text/css" />
  
 
<!------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------->
 
<!------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------->
Line 1,861: Line 1,796:
 
bkgFunction();
 
bkgFunction();
 
</script>
 
</script>
 +
 +
 +
 +
<!----------------------------------------------------------------------------------------------------------------------------------------->
 +
<!---- Left Navigator ---->
 +
<!----------------------------------------------------------------------------------------------------------------------------------------->
 +
 +
 +
<ul class="leftNav" style="margin:0;padding:0;">
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle1">Overview</a>
 +
</li>
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle2">Part I: Yield of recombined P<sub>target</sub></a>
 +
<ul class="leftNavUl2">
 +
<li class="leftNavLi2"><a class="leftNavA2" href="#mainTitle2_1">Induced expression model</a></li>
 +
<li class="leftNavLi2"><a class="leftNavA2" href="#mainTitle2_2">Reverse Transcription model</a></li>
 +
<li class="leftNavLi2"><a class="leftNavA2" href="#mainTitle2_3">Cre Recombination Model</a></li>
 +
<li class="leftNavLi2" style="display:none;"><a class="leftNavA2" href="#mainTitle2_4"></a></li>
 +
</ul>
 +
</li>
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle3">Part II: Times of occurrence of recombined P<sub>target</sub></a>
 +
</li>
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle4">Part III: Optimal induction time</a>
 +
</li>
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle5">Complementary materials</a>
 +
</li>
 +
 +
<li class="leftNavLi"><a class="leftNavA" href="#mainTitle6">Reference</a>
 +
</li>
 +
 +
</ul>
 +
 +
<style>
 +
.leftNav{
 +
position: absolute;
 +
top:20vw;
 +
left:4%;
 +
list-style: none;
 +
z-index: 3;
 +
text-align:left !important;
 +
width:15vw;
 +
}
 +
.leftNav .leftNavA2{
 +
position: relative;
 +
}
 +
.leftNav .leftNavA2:before {
 +
content: "";
 +
display: inline-block;
 +
width: 5px;
 +
height: 5px;
 +
background-color: rgba(55,55,62,0.7);
 +
border-radius: 50%;
 +
position: absolute;
 +
left: 0.3vw;
 +
top:50%;
 +
margin: 0;
 +
padding: 0;
 +
transform: translateY(-50%);
 +
}
 +
.leftNavUl2{
 +
list-style: none;
 +
}
 +
.leftNavA{
 +
display: block;
 +
font-size: 1.1vw;
 +
font-family: Gotham, "Helvetica Neue", Helvetica, Arial, "sans-serif";
 +
padding: 0.2vw 0.8vw;
 +
color:#4CA198;
 +
text-decoration: none;
 +
}
 +
.leftNavLi2{
 +
width: 13vw;
 +
margin:0.3vw  0;
 +
line-height: 17px;
 +
}
 +
.leftNavA:visited{
 +
text-decoration: none;
 +
color:#4CA198;
 +
}
 +
.leftNavA:focus{
 +
text-decoration: none;
 +
color:#4CA198;
 +
}
 +
.leftNavA:hover,.leftNavA:active{
 +
text-decoration: none;
 +
color: #B0EBE4;
 +
}
 +
 +
.leftNavA2{
 +
display: block;
 +
font-size: 1vw;
 +
font-family:Segoe, "Segoe UI", "DejaVu Sans", "Trebuchet MS", Verdana, "sans-serif";
 +
padding: 0.1vw 1.5vw;
 +
color:rgba(168,168,168,1);
 +
text-decoration: none;
 +
}
 +
.leftNavA2:visited{
 +
text-decoration: none;
 +
color:rgba(168,168,168,1);
 +
}
 +
.leftNavA2:focus{
 +
text-decoration: none;
 +
color:rgba(168,168,168,1);
 +
}
 +
.leftNavA2:hover,.leftNavA2:active{
 +
text-decoration: none;
 +
color: white;
 +
}
 +
.menu-active,.menu-active:focus{
 +
text-decoration: none;
 +
color: white;
 +
}
 +
.menu-active2,.menu-active2:focus{
 +
text-decoration:none;
 +
color:#B0EBE4;
 +
}
 +
 +
 +
.mjx-chtml {
 +
    outline: 0;
 +
    font-size:80%;
 +
}
 +
.MJXc-display {
 +
        overflow-x: auto;
 +
        overflow-y: hidden;
 +
}
 +
@media only screen and (max-width:768px){
 +
.leftNav{
 +
display:none;
 +
}
 +
}
 +
</style>
 +
  
 
 
Line 1,867: Line 1,939:
 
<!---- Content ---->
 
<!---- Content ---->
 
<!----------------------------------------------------------------------------------------------------------------------------------------->
 
<!----------------------------------------------------------------------------------------------------------------------------------------->
       
+
            <div class="row title2" id="mainTitle1">
+
<div class="container" id="containerWithLeftNav">
                <div class="col">Modeling</div>
+
<div class="row">
            </div>
+
 +
<div class="row title2" id="mainTitle1">
 +
<div class="col">Overview</div>
 +
</div>
 +
 +
<div class="row para1">
 +
<div class="row">
 +
<div class="col col-lg-12">
 +
Our mutagenesis system uses the BL21 (DE3) <i>E. coli</i> strain transformed with two plasmids, a stringent plasmid named P<sub>target</sub> carrying the target sequence that we want to mutate, and a relaxed plasmid named P<sub>mutant</sub>, carrying the gene encoding the tools necessary for mutagenesis, i.e. reverse transcriptase (RT) and Cre. <br /><br />
 +
As we are designing a brand-new mutagenesis system inside <i>E. coli</i>, we want to demonstrate whether and under what condition it can work, so we turn to modelling to answer these questions. Our modelling work is comprised of 3 parts. 1) We used 3 deterministic models to describe the 3 reaction steps of our system—induced expression, reverse transcription and recombination. This allows us to compute and maximize the yield of the recombined P<sub>target</sub> which in turn, contributes to the optimization of our experimental setup. 2) We simulated the recombination process stochastically and calculated the number of recombined products that occurred during one replication cycle of <i>E. coli</i>. 3) We combined the 3 reaction steps together using deterministic model and found that selecting the least efficient degradation tag for Cre is optimal.
 +
</div>
 +
</div>
 +
</div>
 +
 +
<div class="row title2" id="mainTitle2">
 +
<div class="col">Part I: Deterministic model to compute the yield of recombined P<sub>target</sub></div>
 +
</div>
  
            <div class="row para1">
+
<div class="row para1">
                <div class="col">
+
<div class="col">
                    In our modeling, we successfully simulated the function of our mutagenesis system, and contributed to improve our experimental setup.<br />
+
When we were constructing the plasmid, we encountered a dilemma concerning how RT and Cre should be expressed. Firstly, we thought of putting them both under a same Lac operon so that their expression can be easily induced merely by one kind of inducer—IPTG. Meanwhile, we also considered using different inducers to achieve a more modular design which would be easier to control. As it would take a long time to test which induced expression scheme is better through experiments, we used modelling to test the two constructs. We modelled all the reactions involved and computed the yield of the desired product, i.e. recombined P<sub>target</sub>. Through comparison of the yield acquired using these two induced expression schemes, we decided that the latter scheme should be employed for our system to perform better.  <br /><br />
                    We divided our system into 3 sub-models: induced expression model, reverse transcription model and Cre recombination model. We utilized deterministic and stochastic techniques with parameters derived from our experiments or published papers.<br /><br />
+
By common knowledge we can assume that, if the amount of RT and Cre needs to be different to achieve optimal yield, we should choose the second scheme and put them under different operons. On the contrary, if the yield reaches the maximum under the maximum amount of RT and Cre, the first scheme should be chosen. <br /><br />  
                    Our modeling established theoretical basis for our experiments:<br />
+
In our initial attempt, we found that modelling all the reactions involved is rather difficult, as the reactions are in such a large number and all mixed together. This circumstance makes inspection of the reasonability of our models and parameters impossible. To overcome this issue, we decided to separate these reactions into three minor models and use the steady-state concentration of the substances derived from the previous model as the input of the next model. The three minor models are: <b><i>induced expression model, reverse transcription model and Cre recombination model</i></b>, corresponding to the 3 reaction steps in R-Evolution. The schematic diagram is shown in Fig. 1.
                    1) The determination of simultaneous or separate expression of reverse transcriptase and Cre is guided by models revealing their different working concentration.<br />
+
</div>
                    2) We estimated the optimal expression level and induction time needed to achieve maximal recombination efficiency.<br />
+
</div>
                    3) We modeled the effect of Cre degradation on recombination.<br />
+
                    4) We demonstrated that mutations accumulate accompanying E. coli growth.<br />
+
                    Modeling acted as a shortcut of answering questions concerning experimental setup and revealed new insights into our system. Thus, we believe that our modeling work is very competitive for the best modeling prize.<br /><br />
+
                </div>
+
            </div>
+
  
            <div class="row title2" id="mainTitle1">
+
                                      <div class="row legend">
                <div class="col">Software</div>
+
<div class="row">
            </div>
+
<img src="https://static.igem.org/mediawiki/2019/5/54/T--Fudan-TSI--Fig1.gif" style="width:70%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig1"></a>Figure 1. Workflow of the model.</a></b><br />
 +
<p>Three Grey boxes indicate three major reaction steps in R-Evolution. Arrows indicate the reaction that certain substance is involved. White arrows indicate the case in which substances that originally exist in E.coli act as inputs. Red arrows indicate the case in which intermediates, which are produced in the previous reaction, are generated or involved in next reaction process. The blue arrow indicates the final output that we would like to observe. Inducer – IPTG or aTc (anhydrotetracycline). RT – reverse transcriptase. Cre – Cre recombinase. cDNA – complementary DNA.</p>
 +
</div>
 +
</div>
 +
 +
<div class="row title3" id="mainTitle2_1">
 +
<div class="col"><i>Induced expression model</i></div>
 +
</div>
 +
 +
<div class="row para1">
 +
<div class="col">
 +
We first assumed that both genes encoding RT and Cre are placed together under a lac operon (<a href="#Fig2">Fig 2a</a>). The repressor protein LacI is stably expressed in the cell, 2 molecules of LacI will form a dimer which binds to LacO DNA fragment and represses the expression of RT and Cre. When IPTG is added and transported into the cell, IPTG molecules will bind with LacI and inhibit its binding to LacO. In this way, RT and Cre can be rescued from suppression (<a href="#Ref1">Nikos et al.</a>). The ordinary differential equations (ODEs) describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations can be found in the <a href="#App">appendix</a>.<br /><br />
 +
</div>
 +
</div>
  
            <div class="row para1">
+
<div class="formula para1">
                <div class="col">
+
                        $$
                    Our software simplifies the primer design process for target-specific mutagenesis via reverse transcriptase (RT). We called it tRNA primer designer. Studies have shown that tRNA functions as the primer for in vivo reverse transcription initiation: the 5' end of the tRNA interacts with RT, and the 3' end matches with the mRNA encoding the target.<br />
+
            \begin{aligned}
                    The software consists of 4 parts: reverse transcriptase selection, target sequence input, designed-tRNA visualization, and primer output. Although we only test MMLV RT experimentally, the software can adjust its designing method based on the properties of well-studied RT from 3 species, MMLV, HIV-1 and RSV. Users could design their tRNA primers even for eukaryotic experiments. In addition, we calculate and output the tRNA acceptor stem annealing temperature, as this might be used as an indicator for likelihood to success.
+
            \frac{\text{d}}{\text{d}t}MR &= k_{sMR}\cdot O_{total} - \lambda_{MR}\cdot MR
                </div>
+
            \\
            </div>
+
            \frac{\text{d}}{\text{d}t}R &= k_{sR}\cdot MR-2\cdot k_{2R}\cdot R^2 + 2\cdot k_{2R}R^2- \lambda_R\cdot R
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}R_2 &= k_{2R}\cdot R^2 - k_{-2R}\cdot R_2 - k_r\cdot R_2\cdot O+k_{-r}\cdot (O_{total}-O)-k_{dr1}\cdot R_2\cdot I^2+k_{-dr1}\cdot I_2R_2 -                 \lambda_{R_2}\cdot R_2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}O &= -k_r\cdot R_2\cdot O + k_{-r}\cdot (O_{total}-O) + k_{dr2}\cdot (O_{total}-O)\cdot I^2 - k_{-dr2}\cdot O\cdot I_2R_2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}I &= -2\cdot k_{dr1}\cdot R_2\cdot I^2 + 2\cdot k_{-dr1}\cdot I_2R_2 - 2\cdot k_{dr2}\cdot (O_{total}-O)\cdot I^2 \\\ \ \ \ \ \ \ &\ \ \ + 2\cdot k_{-dr2}\cdot O\cdot I_2R_2 + k_{ft}\cdot YI_{ex} + k_t \cdot (I_{ex}-I) + 2\cdot \lambda_{I_2R_2}\cdot I_2R_2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}I_2R_2 &= k_{dr1}\cdot R_2\cdot I^2 - k_{-dr1} \cdot I_2R_2 + k_{dr2}\cdot (O_{total}-O)\cdot I^2 - k_{-dr2}\cdot O\cdot I_2R_2 - \lambda_{I_2R_2}\cdot I_2R_2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}MY &= k_{sMY}\cdot O_{total} - \lambda_{MY}\cdot MY
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}Y &= k_{sY}\cdot MY + (k_{ft}+k_p)\cdot YI_{ex} - k_p\cdot Y\cdot I_{ex} - \lambda_Y\cdot Y
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}YI_{ex} &= -(k_{ft}+k_p)\cdot YI_{ex} + k_p\cdot Y\cdot I_{ex} - \lambda_{YI_{ex}}\cdot YI_{ex}
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}MRT &= k_{s0RT}\cdot (O_{total}-O) + k_{s1RT}\cdot O - \lambda_{MRT}\cdot MRT
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}RT &= k_{sRT}\cdot MRT - \lambda_{RT}\cdot RT
 +
            \end{aligned}
 +
                        $$
 +
</div>
  
            <div class="row title2" id="mainTitle1">
+
                                        <div class="row legend">
                <div class="col">Hardware</div>
+
<div class="row">
            </div>
+
<img src="https://static.igem.org/mediawiki/2019/6/6f/T--Fudan-TSI--Formula1.gif" style="width:80%; margin:auto;">
 +
</div>
 +
</div>
  
            <div class="row para1">
+
 
                <div class="col">
+
 
                    To track bacteria growth on the plate and observe the fluorescence recovery from nonsense mutation due to continuous mutagenesis, we devised this hardware - the Fluorescence Tracker. It provides continuous, hands-off recording of the growth of plate colonies as well as fluorescent protein expression. For users of our mutagenesis system, with the help of our hardware, they could plate, and then monitor all plates together to increase the likelihood of spotting bacteria colonies with recovered fluorescence at the earliest time point. After discussions with our PI, we improved our hardware by adding remote access through TeamViewer, which allows visualizing the dynamic changes on smartphones. Although the current hardware is only suitable for monitoring fluorescence recovery, it could be easily modified to monitor bacteria colonies growing out of any antibiotic plate. Our hardware allows us to come to lab knowing that a plate with desired colonies is waiting for us.<br />
+
<div class="row para1">
                </div>
+
<div class="col">
            </div>
+
According to our modelling result, the amount of target protein (RT and Cre) will be extremely low when IPTG is not added (<a href="#Fig2">Fig. 2b</a>). The origin point represents the time when an E. coli comes into being through reproduction. As a result, the lac operon is not fully repressed by LacI dimer, causing a leakage expression of target protein (from 0 min to 1 min, <a href="#Fig2">Fig. 2b&c</a>). After that, due to slow degradation rate of the target protein’s mRNA as well as the target protein itself, the amount of target protein will continue to accumulate to a certain amount after the lac operon is fully repressed (from 1 min to 5 min, <a href="#Fig2">Fig. 2b&c</a>). Finally, the degradation process removes target protein from the system (from 5 min to 50 min, <a href="#Fig2">Fig. 2b</a>). When IPTG is added, we find that the concentration of protein product quickly rises as the repression of lac operon is quickly removed (<a href="#Fig2">Fig. 2b&c</a> from 50 min to 100 min ). The steady-state concentration is 6.70 μM. This number will be used for further analysis.
 +
</div>
 +
</div>
 +
 
 +
                                        <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/archive/5/59/20191021145821%21T--Fudan-TSI--Fig2b%26c.gif" style="width:80%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig2">Figure 2. Induced expression of RT and Cre.</a></b><br />
 +
<p><b>a)</b> Schematic diagram of the model. <b>b)</b> Dynamics of target protein. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of protein (RT and Cre) within the system. RT and Cre are expressed under the same Lac operon. <b>c)</b> Dynamics of free lac operon. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of free lac operon, i.e. the lac operon unbound by tetR dimer, within the system. The vertical magenta line indicates the moment when 50μM is added to the system.</p>
 +
</div>
 +
</div>
 +
 
 +
 +
 
 +
<div class="row title3" id="mainTitle2_2">
 +
<div class="col"><i>Reverse Transcription model</i></div>
 +
</div>
 +
 +
<div class="row para1">
 +
<div class="col">
 +
From the first model, the concentration of both RT and Cre are acquired. The concentration of RT serves as input to the reverse transcription model. As the schematic diagram depicts (<a href="#Fig3">Fig. 3a</a>), tRNA primer first binds with reverse transcriptase. When this complex binds with a certain fragment on the target sequence, which is called primer binding site (PBS), the reverse transcription will start and cDNA will be synthesized.<br /><br />
 +
Although a more elaborate model of reverse transcription has been proposed by <a href="#Ref2">Kulpa et al.</a>, it includes many reactions whose kinetic properties are not well characterized. As a result, we simplified that model and came up with our own. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations we used can be found in the <a href="#App">appendix</a>.<br /><br />
 +
</div>
 +
</div>
 +
 
 +
<div class="formula para1">
 +
            $$
 +
            \begin{aligned}
 +
            \frac{\text{d}}{\text{d}t}mGOI &= k_{smGOI}\cdot P_{target} - k_{anneal}\cdot mGOI\cdot C2 - \lambda_{mGOI}\cdot mGOI
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}Pr &= k_{sPr}\cdot P_{mutant} - k_{bind}\cdot RT\cdot Pr + k_{dis}\cdot C2 - \lambda_{Pr}\cdot Pr
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}C2 &= k_{bind}\cdot RT\cdot Pr - \lambda_{C_2}\cdot C2  - k_{anneal}\cdot mGOI\cdot C2 - k_{dis}\cdot C2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}RT &= -k_{bind}\cdot RT+k_{dis}\cdot C2
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}C3 &= k_{anneal}\cdot mGOI\cdot C2 - \lambda_{C3\_RT}\cdot C3
 +
            \\
 +
            \frac{\text{d}}{\text{d}t}cDNA &= k_{scDNA}\cdot C3 - \lambda_{cDNA}\cdot cDNA
 +
            \end{aligned}
 +
            $$
 +
</div>
 +
 
 +
                                        <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/3/3f/T--Fudan-TSI--Formula2.gif" style="width:80%; margin:auto;">
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
The modelling result is shown in <a href="#Fig3">Fig. 3b</a>. It shows that the concentration of cDNA will accumulate at the presence of RT (whose initial concentration is 6.70 μM, computed by the induced expression model) and finally reach a steady-state of 9.60 nM. This number will be used for further analysis.
 +
</div>
 +
</div>
 +
 
 +
                                      <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/3/31/T--Fudan-TSI--Fig3.gif" style="width:80%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig3">Figure 3. Reverse transcription.</a></b> <br />
 +
<p><b>a)</b> Schematic diagram of the model. <b>b)</b> Dynamics of cDNA. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of cDNA within the system.</p>
 +
</div>
 +
</div>
 +
 +
<div class="row title3" id="mainTitle2_3">
 +
<div class="col"><i>Cre Recombination Model</i></div>
 +
</div>
 +
 +
 +
<div class="row para1">
 +
<div class="col">
 +
Our first assumption is that the genes encoding RT and Cre are both placed under lac operon and thus be expressed in the same amount. So now we are about to compute the yield of our desired product to identify whether this experimental setup is feasible. The model of the recombination process has been clearly described by <a href="#Ref3">Ehrilich et al</a>. We made some changes to it according to our own experimental design. The schematic diagram is shown in <a href="#Fig4">Fig. 4a</a>. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations can be found in the <a href="#App">appendix</a>.<br /><br />
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
$$
 +
\begin{aligned}
 +
\frac{\text{d}}{\text{d}t}Ps &= k_{-1}\cdot Ps\_Cre1 - k_1\cdot Ps\cdot Cre - k_{on}\cdot Ps \cdot T7RNA_p + k_{off}\cdot Ps\_T7RNAp
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ds &= k_{-1}\cdot Ds\_Cre1 - k_1\cdot Ds\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ps\_Cre1  &= k_1\cdot Ps\cdot Cre - k_{-1}\cdot Ps\_Cre1 + k_{-2}\cdot Ps\_Cre2 - k_2\cdot Ps\_Cre1\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ds\_Cre1 &= k_1\cdot Ds\cdot Cre - k_{-1}\cdot Ds\_Cre1 + k_{-2}\cdot Ds\_Cre2 - k_2\cdot Ds\_Cre1\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ps\_Cre2 &= -k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 + k_{-34}\cdot Pp\_Dp\_Cre4 + k_2\cdot Ps\_Cre1\cdot Cre - k_{-2} \cdot Ps\_Cre2
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ds\_Cre2 &= -k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 + k_{-34}\cdot Pp\_Dp\_Cre4 + k_2\cdot Ds\_Cre1\cdot Cre - k_{-2}\cdot Ds\_Cre2
 +
\\
 +
\frac{\text{d}}{\text{d}t}Pp\_Dp\_Cre4 &= -k_{-34}\cdot Pp\_Dp\_Cre4 + k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 - k_5\cdot Pp\_Dp\_Cre4 + k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2
 +
\\
 +
\frac{\text{d}}{\text{d}t}Pp\_Cre2 &= k_5\cdot Pp\_Dp\_Cre4 - k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2 + k_2\cdot Pp\_Cre1\cdot Cre - k_{-2}\cdot Pp\_Cre2
 +
\\
 +
\frac{\text{d}}{\text{d}t}Dp\_Cre2 &= k_5\cdot Pp\_Dp\_Cre4 - k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2 + k_2\cdot Dp\_Cre1\cdot Cre - k_{-2}\cdot Dp\_Cre2
 +
\\
 +
\frac{\text{d}}{\text{d}t}Pp\_Cre1 &= - k_2\cdot Pp\_Cre1\cdot Cre + k_{-2}\cdot Pp\_Cre2- k_{-1} \cdot Pp\_Cre1 + k_1\cdot Pp\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Dp\_Cre1 &= - k_2\cdot Dp\_Cre1\cdot Cre + k_{-2}\cdot Dp\_Cre2 - k_{-1} \cdot Dp\_Cre1 + k_1\cdot Dp\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Pp &= k_{-1}\cdot Pp\_Cre1 - k_1\cdot Pp\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Dp &= k_{-1}\cdot Dp\_Cre1 - k_1\cdot Dp\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Cre &= - k_1\cdot Ps\cdot Cre + k_{-1}\cdot Ps\_Cre1 - k_1\cdot Ds\cdot Cre \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ + k_{-1}\cdot Ds\_Cre1 - k_2\cdot Ps\_Cre1\cdot Cre + k_{-2}\cdot Ps\_Cre2 - k_2\cdot Ds\_Cre1\cdot Cre \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ + k_{-2}\cdot Ds\_Cre2 + k_{-2}\cdot Pp\_Cre2 - k_2\cdot Pp\_Cre1\cdot Cre + k_{-2}\cdot Dp\_Cre2 \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ - k_2\cdot Dp\_Cre1\cdot Cre + k_{-1}\cdot Pp\_Cre1 - k_1\cdot Pp\cdot Cre + k_{-1}\cdot Dp\_Cre1 - k_1\cdot Dp\cdot Cre
 +
\\
 +
\frac{\text{d}}{\text{d}t}Ps\_T7RNAp &= k_{on}\cdot Ps\cdot T7RNAp - k_{off}\cdot Ps\_T7RNAp
 +
\\
 +
\frac{\text{d}}{\text{d}t}T7RNAp &= 0
 +
\end{aligned}
 +
$$
 +
</div>
 +
 
 +
                                        <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/e/e9/T--Fudan-TSI--Formula3.gif" style="width:80%; margin:auto;">
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
As is shown in the diagram, 2 Cre molecules bind with 1 loxP site successively, either on cDNA or P<sub>target</sub>. Four Cre molecules will form a Holliday junction, and thus starting the recombination reaction. Two pairs of loxP will work together and complete the strand exchange between cDNA and P<sub>target</sub>. After that, the recombined product is produced. What we are interested in is the percentage of recombined P<sub>target</sub> among all P<sub>target</sub>s in one E. coli. So, we turn to compute that percentage based on the model that we have established.<br /><br />
 +
Unfortunately, we found that the amount of substances is too small. For example, the concentration of  is only 10 nM, which means there are only about 5 molecules of  in one cell. These small numbers caused some computational problems in Matlab when we were using its ODE solver (ode15s). To address this problem, we converted the units of the amount of the substances from mole per litter (M) to molecule. The units of the kinetic parameters are also converted accordingly.<br /><br />
 +
Now the recombination step is modeled under the initial condition of 5 molecules of non-mutated , 3228 molecules of Cre and 5 molecules of cDNA (<a href="#Fig4">Fig 4b</a>). The last two numbers are the outputs of previous models after going through some unit conversion steps.
 +
</div>
 +
</div>
 +
<div class="row para1">
 +
<div class="col">
 +
The result is disappointing. After a long period of reaction, no recombined  showed up. It is because there are too many Cre molecules so s are all bounded by them and remain in the intermediate form. What’s more,  can't bind with T7 RNA polymerase and be transcribed as a consequence of Cre occupation. This leads to the system’s inability of undergoing further reverse transcription process, stopping cDNA’s production, resulting in a stop of the system, and rendering mutation accumulation impossible (<a href="#Fig4">Fig. 4c</a>).<br /><br />
 +
                                                      This result tells us that the number of Cre molecules needs to be much lower for the system to function. We then set out to determine how many Cre is optimal. After we fed the recombination model with cDNA and Cre at different concentrations, the problem seems to be clear as the yield of recombined  varies greatly responding to different numbers of cDNA and Cre (<a href="#Fig4">Fig. 4d</a>). When cDNA is confined to 5 molecules, we will get no yield at all in the period of E. coli's replication cycle if the concentration of Cre is greater than 20 molecules. Instead, the yield is maximized when the final Cre concentration is around 2 molecules (<a href="#Fig4">Fig 4e</a>).
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
Now we use the optimized number of Cre as the input to our third model. The result is shown in <a href="#Fig4">Fig. 4f</a>, which is satisfactory. The recombined P<sub>target</sub>) finally occurs and  has a chance to bind with T7 RNA polymerase, which means mutated gene of interest could be transcribed and further mutated, thus making the accumulation of mutations possible (Fig 4g). These results remind us to use different inducer to induce the expression of RT and Cre. So, we revised our experimental design and decided to use Tet operon to control the expression of Cre and induce that with anhydrotetracycline (aTc). Even though we later used degradation tag to accelerate the degradation process of Cre and to decrease the expression level of Cre, considering the fact that the tet operon is less prone to leakage and that using merely lac operon to control the expression of RT and Cre may cause unexpected problems, we still used different operons to control the expression of RT and Cre. This setup will be considered in the model in Part III.
 +
</div>
 +
</div>
 +
<div class="row para1">
 +
<div class="col">
 +
There is still something that is not well explained in our current model. The final percentage of recombined P<sub>target</sub> is around 1.5%. The unit of the substance is molecules, so it means there is 0.075 recombined P<sub>target</sub> in one cell, which is unrealistic. This problem reflects that converting the unit of substance into molecule when doing deterministic modelling cannot offer a precise description of the system’s status.<br /><br />
 +
We then used stochastic modelling techniques to determine whether and how many recombined P<sub>target</sub>s will show up in one replication cycle of E. coli.
 +
</div>
 +
</div>
 +
 +
                                      <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/d/d8/T--Fudan-TSI--Fig4.gif" style="width:70%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig4">Figure 4. Cre recombination (deterministic).</a></b><br />
 +
<p><b>a)</b> Schematic diagram of the model. Ps, Un-recombined P<sub>target</sub>. Pp, Recombined P<sub>target</sub>. <b>b-c)</b> Recombination when Cre is expressed under Lac operon. Dynamics of the percentage of un-recombined/ recombined P<sub>target</sub> among all P<sub>target</sub>s is shown in <b>b</b>. Horizontal axis shows the length of time (8 hours, corresponding to R-Evolution’s function period). The distribution of the percentage of substances at the steady-state is shown in <b>c</b>. <b>d)</b>  Yield of recombined P<sub>target</sub> at different initial number of cDNA and Cre. The yield of recombined P<sub>target</sub> is calculated as the percentage of recombined P<sub>target</sub> among all P<sub>target</sub>s. The horizontal white line corresponds to current situation where the initial number of cDNA is 5 molecules in one E.coli. <b>e)</b>  Yield of recombined P<sub>target</sub> at different initial number of Cre when initial number of cDNA is 5 molecules. <b>f-g)</b>  Recombination when Cre is expressed under different operon. Dynamics of the percentage of un-recombined/ recombined P<sub>target</sub> among all P<sub>target</sub>s is shown in <b>f</b>. The distribution of the percentage of substances at the steady-state is shown in <b>g</b>.</p>
 +
 
 +
</div>
 +
</div>
 +
 
 +
<div class="row title3" id="mainTitle2_4">
 +
</div>
 +
 +
<div class="row title2" id="mainTitle3">
 +
<div class="col">Part II: Stochastic model to compute times of occurrence of recombined P<sub>target</sub></div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
We use Gillespie algorithm in stochastic modelling. The procedure of this algorithm is shown as follows in the form of pseudocode:<br /><br />
 +
                            Step 1: Initialize the reaction system at \(t=0\) with rate constants \(c_1, c_2, ......, c_v\) as initial numbers of molecules \(x_1, x_2, ......, x_u\) corresponding to \(v\) reactions and \(u\) sustances (both reactants and products) involved in the reaction system.<br /><br />
 +
                            Step 2: For each \(i=1,2,......,v\), calculate the hazard for the \(i^{th}\) type of reaction, denoted as \(h_i(x,c_i)\) based on current substance number \(x\).<br /><br />
 +
                            Step 3: Calculate the combined reaction hazard \(h_0(x,c)=\Sigma_{i=1}^{v}h_i(x,c_i)\).<br /><br />
 +
                            Step 4: Simulate the time to the next reaction, \(t^{'}\) , which is a random quantity that obeys exponential distribution with parameter \(\lambda\).<br /><br />
 +
                            Step 5: Put \(t:=t+t^{'}\). <br /><br />
 +
                            Step 6: Simulate the reaction index, \(j\). The probability that the \(j^{th}\) reaction can occur is \(\frac{h_i(x,c_i)}{h_0(x,c)}, i=1,2,......,v.\).<br /><br />
 +
                            Step 7: Update \(x\) according to reaction \(j\), which means putting \(x:=x+S^{(j)}\), where \(S^{(j)}\) denotes the \(j^{th}\) colomn of the stoichiometry matrix \(S\). The \(j^{th}\) column of  denotes the change in substance number  caused by the \(j^{th}\) reaction.<br /><br />
 +
                            Step 8: Record time \(t\) and current substance number \(x\).<br /><br />
 +
                            Step 9: If \(t\ \text{<}\ T_{max}\), return to step 2. \(T_{max}\) corresponds to the maximum duration of the reaction set by the user.<br /><br />
 +
                            Step 10: Plot the result to see the dynamic of the quantity of the substance that you are interested in.<br /><br />
 +
                            Although the algorithm is rather simple, basic mathematical skills is required to understand its theoretical basis. You may consult the book written by <a href="#Ref4">Wilkinson and Darren J.</a> for a thorough understanding. The result is shown in <a href="#Fig5">Fig. 5</a>.<br /><br />
 +
                            The result demonstrates that recombined P<sub>target</sub>s do occur and two rounds of reverse transcription and recombination can take place in one replication cycle of E. coli (1200 s) (<a href="#Fig5">Fig 5a</a>). On the contrary, no recombined  will come out within that period if the initial cDNA is 5 molecules and initial Cre is 3228 molecules (<a href="#Fig5">Fig 5b</a>). This again demonstrates the necessity of putting RT and Cre under different induction setups. The fluctuation of the number of recombined P<sub>target</sub> is due to the backward reaction that Cre can rebind with recombined  and reverting the action, making it not counted as recombined P<sub>target</sub> by the algorithm.
 +
</div>
 +
</div>
 +
 
 +
                                      <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/5/5f/T--Fudan-TSI--Fig5bc.gif" style="width:80%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig5">Figure 5. Cre recombination (stochastic).</a></b><br />
 +
<p>Horizontal axis shows the length of time (20min, corresponding to the duration of 1 E.coli replication cycle). Vertical axis demonstrates the number of recombined P<sub>target</sub>. The initial number of Cre is 2 molecules in <b>a</b>, 3228 molecules in <b>b</b>.</p>
 +
</div>
 +
</div>
 +
 +
<div class="row title2" id="mainTitle4">
 +
<div class="col">Part III: Deterministic model to determine optimal induction time</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
        To ensure the evolved protein is encoded by the mutated GOI sequence that is recombined into P<sub>target</sub>, we decided to use degradation tag to accelerate the degradation process of Cre. This design would make Cre only function when inducer is in the system, thus allowing stringent control of the protein. However, we then face the problem of how to select the optimal degradation tag. Empirically, to minimize the duration of recombination, we tend to choose degradation tags with higher efficiency, but extremely high degradation rate will also reduce the yield of recombined P<sub>target</sub>, leading to decreased library size. Also, it is impractical for researchers to do experiments to test these degradation tags one by one. For these reasons, we are going to use models to find out the optimal degradation tag that should be added to Cre based on the average yield of recombined P<sub>target</sub> at the end of R-Evolution functioning period (8 hours).<br /><br />
 +
                            We intend to use the models described in Part I, combined with aTc induction model proposed by <a href="#Ref5">Steel et al</a>. to compute the yield of recombined P<sub>target</sub> under different degradation rate of Cre (the reason why Tet operon is used has been elaborated in Part I; the schematic diagram of this process is shown in <a href="#Fig6">Fig 6a</a>. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations we used can be found in the <a href="#App">appendix</a>.
 +
</div>
 +
</div>
 +
 
 +
                                      <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/9/9c/T--Fudan-TSI--Formula4.gif" style="width:80%; margin:auto;">
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
                            Although the setup in Part I successfully provided us with a clear insight into the reactions and dynamic changes of substances that underlie our mutagenesis system, the simplification that the steady-state substance concentrations of previous models can be used as inputs for subsequent models doesn’t match real reaction situation. For example, when Cre is expressed, it can immediately bind with cDNA and initiate recombination. This fact contradicts with our model assumption that recombination only takes place after both Cre and cDNA has reached their steady-state concentration.<br /><br />
 +
                            To overcome this issue, we decided to combine all three minor models together and calculate the expected output.<br /><br />
 +
                            As a result of the impreciseness of the basic assumption of the models in part I, we only gave a qualitative conclusion that the amount of RT and Cre should be different. Here we need to quantify how Cre degradation rate and steady-state concentration affects the yield of recombined P<sub>target</sub>. That’s why we employed deterministic model here to combine the separate steps together into one and better simulate real intracellular circumstances.<br /><br />
 +
                            By combining the models that have been talked above, we revealed the reason why the degradation tag with a moderate degradation rate, which can’t be too high or too low, should be selected (<a href="#Fig6">Fig 6b</a>): under appropriate inducer concentration (20~22uM), when the degradation rate is relatively low (below 0.1 min<sup>-1</sup>), the yield of recombined P<sub>target</sub> will increase according to the increase of Cre degradation rate, but when that rate is sufficiently high (above 0.1 min<sup>-1</sup>), the increase of Cre degradation rate will do harm to the yield of recombined P<sub>target</sub>.<br /><br />
 +
                            The average degradation rate acquired from literature is 0.2 min<sup>-1</sup>(<a href="#Ref1">Nikos et al.</a>) and the degradation rate of Cre when tagged with the most efficient degradation tag is 0.69 min<sup>-1</sup>. Within this range of degradation rate, the maximum yield of recombined P<sub>target</sub> will decrease according to the increase of Cre degradation efficiency (<a href="#Fig6">Fig 6c</a>). So we decided to use the least efficient degradation tag.<br /><br />
 +
                            We also revealed the dynamic change of the recombined P<sub>target</sub>. It will continuously accumulate within Cre function period (<a href="#Fig6">Fig 6d</a>). However, the concentration remains to be low within that period, due to Cre degradation (<a href="#Fig6">Fig 6e</a>).<br /><br />
 +
                            Finally, there is another interesting phenomenon that is worth mentioning. From <a href="#Fig6">Fig 6b</a> and <a href="#Fig6">Fig 6c</a>, we can find that for each degradation tag rate greater than 0.2 min<sup>-1</sup>, there exits a range of aTc dosage that can make the yield of recombined  relatively big. Also, decreased degradation efficiency enlarges that range. This discovery provides us with another reason for using less efficient degradation tag in that it can increase the robustness of our mutagenesis system by decreasing its sensitivity to the change of inducer dosage.
 +
</div>
 +
</div>
 +
 
 +
                                      <div class="row legend">
 +
<div class="row">
 +
<img src="https://static.igem.org/mediawiki/2019/archive/8/84/20191020042934%21T--Fudan-TSI--Fig6.gif" style="width:80%; margin:auto;">
 +
</div>
 +
<div class="row legends">
 +
<b><a name="Fig6">Figure 6. Whole process simulation considering Cre degradation tag.</a></b> <br />  
 +
<p><b>a)</b> Chemical reactions involved in Cre induced expression. <b>b)</b> Yield of recombined P<sub>target</sub> at different Cre degradation rate and aTc dosage. The white line on the left corresponding to the case where the degradation rate of Cre is 0.2 min<sup>-1</sup>. The white line on the left corresponding to the case where the degradation rate of Cre is 0.69 min<sup>-1</sup>. <b>c)</b> Yield of recombined P<sub>target</sub> at different Cre degradation rate and aTc dosage (3D plot). The range of Cre degradation rate is 0.2~0.69 min<sup>-1</sup>. <b>d)</b> Dynamics of yield of recombined P<sub>target</sub> at the Cre degradation rate of 0.2 min<sup>-1</sup> and the initial 22uM aTc dosage. <b>e)</b> Dynamics of yield of recombined P<sub>target</sub> at the Cre degradation rate of 0.2 min<sup>-1</sup> and the initial 22uM aTc dosage.</p>
 +
</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
In the deterministic model, we combined the three minor models proposed previously and assessed the mutagenesis system in whole. Through this addition, we achieved a better simulation of the real intracellular reactions and answered the question of when Cre should be induced for the highest level of recombination efficiency to be obtained.
 +
</div>
 +
</div>
 +
 +
<div class="row title3" id="mainTitle5">
 +
<div class="col">Complementary materials</div>
 +
</div>
 +
 
 +
<div class="row para1">
 +
<div class="col">
 +
You can access the Matlab codes of our models from our <a href="https://github.com/kiyochou/2019-iGEM-Fudan-TSI" target="_blank">Github</a> repository.<br /><br />
 +
Please consult the appendix file for a clearer understanding of the formulation of the model.<br /><br />
 +
<a name="App"><embed src="https://static.igem.org/mediawiki/2019/8/8c/T--Fudan-TSI--ModelApp.pdf" width="100%" height="1000"  type="application/pdf"></a>
 +
</div>
 +
</div>
 +
 
 +
<div class="row title3" id="mainTitle6">
 +
<div class="col">References</div>
 +
</div>
 +
 
 +
 +
<div class="row para1">
 +
<div class="col">
 +
<ul class="paraUl" style="list-style:none;">
 +
<li>[1]. <a name="Ref1" href="https://www.cell.com/biophysj/fulltext/S0006-3495(08)00100-8" target="_blank">Stamatakis M, Mantzaris N V. Comparison of Deterministic and Stochastic Models of the lac Operon Genetic Network[J]. Biophysical Journal, 2009, 96(3):887-906.</a></li>
 +
<li>[2]. <a name="Ref2" href="https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1169686/pdf/000856.pdf" target="_blank">Kulpa, D. Determination of the site of first strand transfer during Moloney murine leukemia virus reverse transcription and identification of strand transfer-associated reverse transcriptase errors[J]. EMBO (European Molecular Biology Organization) Journal, 1997, 16(4):856-865.</a></li>
 +
<li>[3]. <a name="Ref3" href="https://www.sciencedirect.com/science/article/pii/S0022283698921490" target="_blank">Ringrose L, Lounnas V, Ehrlich L, et al. Comparative kinetic analysis of FLP and cre recombinases: mathematical models for DNA binding and recombination[J]. Journal of Molecular Biology, 1998, 284(2):0-384.</a></li>
 +
<li>[4]. <a name="Ref4">Wilkinson, Darren J. Stochastic Modelling for Systems Biology, Second Edition[M]. Crc Press, 2011.</a></li>
 +
<li>[5]. <a name="Ref5" href="https://ieeexplore.ieee.org/document/8263882" target="_blank">Harris A W K, Kelly C L, Steel H, et al. The autorepressor: A case study of the importance of model selection[C]. Decision & Control. IEEE, 2018.</a></li>
 +
</ul>
 +
</div>
 +
</div>
 +
 
 +
 +
</div>
 +
</div>
 
 
  
Line 2,243: Line 2,635:
 
<!---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------->
 
<!---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------->
  
 
</p>
 
<!--
 
NewPP limit report
 
CPU time usage: 0.010 seconds
 
Real time usage: 0.011 seconds
 
Preprocessor visited node count: 6/1000000
 
Preprocessor generated node count: 34/1000000
 
Post‐expand include size: 0/2097152 bytes
 
Template argument size: 0/2097152 bytes
 
Highest expansion depth: 2/40
 
Expensive parser function count: 0/100
 
-->
 
 
<!-- Saved in parser cache with key 2019_igem_org:pcache:idhash:6656-0!*!*!*!*!*!* and timestamp 20191021203154 and revision id 489243
 
-->
 
</div>             <div class="visualClear"></div>
 
            </div>
 
    </div>
 
 
        </div>
 
    </div>
 
 
</html>
 
</html>

Revision as of 21:37, 21 October 2019

Overview
Our mutagenesis system uses the BL21 (DE3) E. coli strain transformed with two plasmids, a stringent plasmid named Ptarget carrying the target sequence that we want to mutate, and a relaxed plasmid named Pmutant, carrying the gene encoding the tools necessary for mutagenesis, i.e. reverse transcriptase (RT) and Cre.

As we are designing a brand-new mutagenesis system inside E. coli, we want to demonstrate whether and under what condition it can work, so we turn to modelling to answer these questions. Our modelling work is comprised of 3 parts. 1) We used 3 deterministic models to describe the 3 reaction steps of our system—induced expression, reverse transcription and recombination. This allows us to compute and maximize the yield of the recombined Ptarget which in turn, contributes to the optimization of our experimental setup. 2) We simulated the recombination process stochastically and calculated the number of recombined products that occurred during one replication cycle of E. coli. 3) We combined the 3 reaction steps together using deterministic model and found that selecting the least efficient degradation tag for Cre is optimal.
Part I: Deterministic model to compute the yield of recombined Ptarget
When we were constructing the plasmid, we encountered a dilemma concerning how RT and Cre should be expressed. Firstly, we thought of putting them both under a same Lac operon so that their expression can be easily induced merely by one kind of inducer—IPTG. Meanwhile, we also considered using different inducers to achieve a more modular design which would be easier to control. As it would take a long time to test which induced expression scheme is better through experiments, we used modelling to test the two constructs. We modelled all the reactions involved and computed the yield of the desired product, i.e. recombined Ptarget. Through comparison of the yield acquired using these two induced expression schemes, we decided that the latter scheme should be employed for our system to perform better.

By common knowledge we can assume that, if the amount of RT and Cre needs to be different to achieve optimal yield, we should choose the second scheme and put them under different operons. On the contrary, if the yield reaches the maximum under the maximum amount of RT and Cre, the first scheme should be chosen.

In our initial attempt, we found that modelling all the reactions involved is rather difficult, as the reactions are in such a large number and all mixed together. This circumstance makes inspection of the reasonability of our models and parameters impossible. To overcome this issue, we decided to separate these reactions into three minor models and use the steady-state concentration of the substances derived from the previous model as the input of the next model. The three minor models are: induced expression model, reverse transcription model and Cre recombination model, corresponding to the 3 reaction steps in R-Evolution. The schematic diagram is shown in Fig. 1.
Figure 1. Workflow of the model.

Three Grey boxes indicate three major reaction steps in R-Evolution. Arrows indicate the reaction that certain substance is involved. White arrows indicate the case in which substances that originally exist in E.coli act as inputs. Red arrows indicate the case in which intermediates, which are produced in the previous reaction, are generated or involved in next reaction process. The blue arrow indicates the final output that we would like to observe. Inducer – IPTG or aTc (anhydrotetracycline). RT – reverse transcriptase. Cre – Cre recombinase. cDNA – complementary DNA.

Induced expression model
We first assumed that both genes encoding RT and Cre are placed together under a lac operon (Fig 2a). The repressor protein LacI is stably expressed in the cell, 2 molecules of LacI will form a dimer which binds to LacO DNA fragment and represses the expression of RT and Cre. When IPTG is added and transported into the cell, IPTG molecules will bind with LacI and inhibit its binding to LacO. In this way, RT and Cre can be rescued from suppression (Nikos et al.). The ordinary differential equations (ODEs) describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations can be found in the appendix.

$$ \begin{aligned} \frac{\text{d}}{\text{d}t}MR &= k_{sMR}\cdot O_{total} - \lambda_{MR}\cdot MR \\ \frac{\text{d}}{\text{d}t}R &= k_{sR}\cdot MR-2\cdot k_{2R}\cdot R^2 + 2\cdot k_{2R}R^2- \lambda_R\cdot R \\ \frac{\text{d}}{\text{d}t}R_2 &= k_{2R}\cdot R^2 - k_{-2R}\cdot R_2 - k_r\cdot R_2\cdot O+k_{-r}\cdot (O_{total}-O)-k_{dr1}\cdot R_2\cdot I^2+k_{-dr1}\cdot I_2R_2 - \lambda_{R_2}\cdot R_2 \\ \frac{\text{d}}{\text{d}t}O &= -k_r\cdot R_2\cdot O + k_{-r}\cdot (O_{total}-O) + k_{dr2}\cdot (O_{total}-O)\cdot I^2 - k_{-dr2}\cdot O\cdot I_2R_2 \\ \frac{\text{d}}{\text{d}t}I &= -2\cdot k_{dr1}\cdot R_2\cdot I^2 + 2\cdot k_{-dr1}\cdot I_2R_2 - 2\cdot k_{dr2}\cdot (O_{total}-O)\cdot I^2 \\\ \ \ \ \ \ \ &\ \ \ + 2\cdot k_{-dr2}\cdot O\cdot I_2R_2 + k_{ft}\cdot YI_{ex} + k_t \cdot (I_{ex}-I) + 2\cdot \lambda_{I_2R_2}\cdot I_2R_2 \\ \frac{\text{d}}{\text{d}t}I_2R_2 &= k_{dr1}\cdot R_2\cdot I^2 - k_{-dr1} \cdot I_2R_2 + k_{dr2}\cdot (O_{total}-O)\cdot I^2 - k_{-dr2}\cdot O\cdot I_2R_2 - \lambda_{I_2R_2}\cdot I_2R_2 \\ \frac{\text{d}}{\text{d}t}MY &= k_{sMY}\cdot O_{total} - \lambda_{MY}\cdot MY \\ \frac{\text{d}}{\text{d}t}Y &= k_{sY}\cdot MY + (k_{ft}+k_p)\cdot YI_{ex} - k_p\cdot Y\cdot I_{ex} - \lambda_Y\cdot Y \\ \frac{\text{d}}{\text{d}t}YI_{ex} &= -(k_{ft}+k_p)\cdot YI_{ex} + k_p\cdot Y\cdot I_{ex} - \lambda_{YI_{ex}}\cdot YI_{ex} \\ \frac{\text{d}}{\text{d}t}MRT &= k_{s0RT}\cdot (O_{total}-O) + k_{s1RT}\cdot O - \lambda_{MRT}\cdot MRT \\ \frac{\text{d}}{\text{d}t}RT &= k_{sRT}\cdot MRT - \lambda_{RT}\cdot RT \end{aligned} $$
According to our modelling result, the amount of target protein (RT and Cre) will be extremely low when IPTG is not added (Fig. 2b). The origin point represents the time when an E. coli comes into being through reproduction. As a result, the lac operon is not fully repressed by LacI dimer, causing a leakage expression of target protein (from 0 min to 1 min, Fig. 2b&c). After that, due to slow degradation rate of the target protein’s mRNA as well as the target protein itself, the amount of target protein will continue to accumulate to a certain amount after the lac operon is fully repressed (from 1 min to 5 min, Fig. 2b&c). Finally, the degradation process removes target protein from the system (from 5 min to 50 min, Fig. 2b). When IPTG is added, we find that the concentration of protein product quickly rises as the repression of lac operon is quickly removed (Fig. 2b&c from 50 min to 100 min ). The steady-state concentration is 6.70 μM. This number will be used for further analysis.
Figure 2. Induced expression of RT and Cre.

a) Schematic diagram of the model. b) Dynamics of target protein. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of protein (RT and Cre) within the system. RT and Cre are expressed under the same Lac operon. c) Dynamics of free lac operon. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of free lac operon, i.e. the lac operon unbound by tetR dimer, within the system. The vertical magenta line indicates the moment when 50μM is added to the system.

Reverse Transcription model
From the first model, the concentration of both RT and Cre are acquired. The concentration of RT serves as input to the reverse transcription model. As the schematic diagram depicts (Fig. 3a), tRNA primer first binds with reverse transcriptase. When this complex binds with a certain fragment on the target sequence, which is called primer binding site (PBS), the reverse transcription will start and cDNA will be synthesized.

Although a more elaborate model of reverse transcription has been proposed by Kulpa et al., it includes many reactions whose kinetic properties are not well characterized. As a result, we simplified that model and came up with our own. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations we used can be found in the appendix.

$$ \begin{aligned} \frac{\text{d}}{\text{d}t}mGOI &= k_{smGOI}\cdot P_{target} - k_{anneal}\cdot mGOI\cdot C2 - \lambda_{mGOI}\cdot mGOI \\ \frac{\text{d}}{\text{d}t}Pr &= k_{sPr}\cdot P_{mutant} - k_{bind}\cdot RT\cdot Pr + k_{dis}\cdot C2 - \lambda_{Pr}\cdot Pr \\ \frac{\text{d}}{\text{d}t}C2 &= k_{bind}\cdot RT\cdot Pr - \lambda_{C_2}\cdot C2 - k_{anneal}\cdot mGOI\cdot C2 - k_{dis}\cdot C2 \\ \frac{\text{d}}{\text{d}t}RT &= -k_{bind}\cdot RT+k_{dis}\cdot C2 \\ \frac{\text{d}}{\text{d}t}C3 &= k_{anneal}\cdot mGOI\cdot C2 - \lambda_{C3\_RT}\cdot C3 \\ \frac{\text{d}}{\text{d}t}cDNA &= k_{scDNA}\cdot C3 - \lambda_{cDNA}\cdot cDNA \end{aligned} $$
The modelling result is shown in Fig. 3b. It shows that the concentration of cDNA will accumulate at the presence of RT (whose initial concentration is 6.70 μM, computed by the induced expression model) and finally reach a steady-state of 9.60 nM. This number will be used for further analysis.
Figure 3. Reverse transcription.

a) Schematic diagram of the model. b) Dynamics of cDNA. Horizontal axis shows the length of time. Vertical axis demonstrates the amount of cDNA within the system.

Cre Recombination Model
Our first assumption is that the genes encoding RT and Cre are both placed under lac operon and thus be expressed in the same amount. So now we are about to compute the yield of our desired product to identify whether this experimental setup is feasible. The model of the recombination process has been clearly described by Ehrilich et al. We made some changes to it according to our own experimental design. The schematic diagram is shown in Fig. 4a. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations can be found in the appendix.

$$ \begin{aligned} \frac{\text{d}}{\text{d}t}Ps &= k_{-1}\cdot Ps\_Cre1 - k_1\cdot Ps\cdot Cre - k_{on}\cdot Ps \cdot T7RNA_p + k_{off}\cdot Ps\_T7RNAp \\ \frac{\text{d}}{\text{d}t}Ds &= k_{-1}\cdot Ds\_Cre1 - k_1\cdot Ds\cdot Cre \\ \frac{\text{d}}{\text{d}t}Ps\_Cre1 &= k_1\cdot Ps\cdot Cre - k_{-1}\cdot Ps\_Cre1 + k_{-2}\cdot Ps\_Cre2 - k_2\cdot Ps\_Cre1\cdot Cre \\ \frac{\text{d}}{\text{d}t}Ds\_Cre1 &= k_1\cdot Ds\cdot Cre - k_{-1}\cdot Ds\_Cre1 + k_{-2}\cdot Ds\_Cre2 - k_2\cdot Ds\_Cre1\cdot Cre \\ \frac{\text{d}}{\text{d}t}Ps\_Cre2 &= -k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 + k_{-34}\cdot Pp\_Dp\_Cre4 + k_2\cdot Ps\_Cre1\cdot Cre - k_{-2} \cdot Ps\_Cre2 \\ \frac{\text{d}}{\text{d}t}Ds\_Cre2 &= -k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 + k_{-34}\cdot Pp\_Dp\_Cre4 + k_2\cdot Ds\_Cre1\cdot Cre - k_{-2}\cdot Ds\_Cre2 \\ \frac{\text{d}}{\text{d}t}Pp\_Dp\_Cre4 &= -k_{-34}\cdot Pp\_Dp\_Cre4 + k_{34}\cdot Ps\_Cre2\cdot Ds\_Cre2 - k_5\cdot Pp\_Dp\_Cre4 + k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2 \\ \frac{\text{d}}{\text{d}t}Pp\_Cre2 &= k_5\cdot Pp\_Dp\_Cre4 - k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2 + k_2\cdot Pp\_Cre1\cdot Cre - k_{-2}\cdot Pp\_Cre2 \\ \frac{\text{d}}{\text{d}t}Dp\_Cre2 &= k_5\cdot Pp\_Dp\_Cre4 - k_{-5}\cdot Pp\_Cre2\cdot Dp\_Cre2 + k_2\cdot Dp\_Cre1\cdot Cre - k_{-2}\cdot Dp\_Cre2 \\ \frac{\text{d}}{\text{d}t}Pp\_Cre1 &= - k_2\cdot Pp\_Cre1\cdot Cre + k_{-2}\cdot Pp\_Cre2- k_{-1} \cdot Pp\_Cre1 + k_1\cdot Pp\cdot Cre \\ \frac{\text{d}}{\text{d}t}Dp\_Cre1 &= - k_2\cdot Dp\_Cre1\cdot Cre + k_{-2}\cdot Dp\_Cre2 - k_{-1} \cdot Dp\_Cre1 + k_1\cdot Dp\cdot Cre \\ \frac{\text{d}}{\text{d}t}Pp &= k_{-1}\cdot Pp\_Cre1 - k_1\cdot Pp\cdot Cre \\ \frac{\text{d}}{\text{d}t}Dp &= k_{-1}\cdot Dp\_Cre1 - k_1\cdot Dp\cdot Cre \\ \frac{\text{d}}{\text{d}t}Cre &= - k_1\cdot Ps\cdot Cre + k_{-1}\cdot Ps\_Cre1 - k_1\cdot Ds\cdot Cre \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ + k_{-1}\cdot Ds\_Cre1 - k_2\cdot Ps\_Cre1\cdot Cre + k_{-2}\cdot Ps\_Cre2 - k_2\cdot Ds\_Cre1\cdot Cre \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ + k_{-2}\cdot Ds\_Cre2 + k_{-2}\cdot Pp\_Cre2 - k_2\cdot Pp\_Cre1\cdot Cre + k_{-2}\cdot Dp\_Cre2 \\\ \ \ \ \ \ \ \ \ \ \ \ &\ \ \ - k_2\cdot Dp\_Cre1\cdot Cre + k_{-1}\cdot Pp\_Cre1 - k_1\cdot Pp\cdot Cre + k_{-1}\cdot Dp\_Cre1 - k_1\cdot Dp\cdot Cre \\ \frac{\text{d}}{\text{d}t}Ps\_T7RNAp &= k_{on}\cdot Ps\cdot T7RNAp - k_{off}\cdot Ps\_T7RNAp \\ \frac{\text{d}}{\text{d}t}T7RNAp &= 0 \end{aligned} $$
As is shown in the diagram, 2 Cre molecules bind with 1 loxP site successively, either on cDNA or Ptarget. Four Cre molecules will form a Holliday junction, and thus starting the recombination reaction. Two pairs of loxP will work together and complete the strand exchange between cDNA and Ptarget. After that, the recombined product is produced. What we are interested in is the percentage of recombined Ptarget among all Ptargets in one E. coli. So, we turn to compute that percentage based on the model that we have established.

Unfortunately, we found that the amount of substances is too small. For example, the concentration of is only 10 nM, which means there are only about 5 molecules of in one cell. These small numbers caused some computational problems in Matlab when we were using its ODE solver (ode15s). To address this problem, we converted the units of the amount of the substances from mole per litter (M) to molecule. The units of the kinetic parameters are also converted accordingly.

Now the recombination step is modeled under the initial condition of 5 molecules of non-mutated , 3228 molecules of Cre and 5 molecules of cDNA (Fig 4b). The last two numbers are the outputs of previous models after going through some unit conversion steps.
The result is disappointing. After a long period of reaction, no recombined showed up. It is because there are too many Cre molecules so s are all bounded by them and remain in the intermediate form. What’s more, can't bind with T7 RNA polymerase and be transcribed as a consequence of Cre occupation. This leads to the system’s inability of undergoing further reverse transcription process, stopping cDNA’s production, resulting in a stop of the system, and rendering mutation accumulation impossible (Fig. 4c).

This result tells us that the number of Cre molecules needs to be much lower for the system to function. We then set out to determine how many Cre is optimal. After we fed the recombination model with cDNA and Cre at different concentrations, the problem seems to be clear as the yield of recombined varies greatly responding to different numbers of cDNA and Cre (Fig. 4d). When cDNA is confined to 5 molecules, we will get no yield at all in the period of E. coli's replication cycle if the concentration of Cre is greater than 20 molecules. Instead, the yield is maximized when the final Cre concentration is around 2 molecules (Fig 4e).
Now we use the optimized number of Cre as the input to our third model. The result is shown in Fig. 4f, which is satisfactory. The recombined Ptarget) finally occurs and has a chance to bind with T7 RNA polymerase, which means mutated gene of interest could be transcribed and further mutated, thus making the accumulation of mutations possible (Fig 4g). These results remind us to use different inducer to induce the expression of RT and Cre. So, we revised our experimental design and decided to use Tet operon to control the expression of Cre and induce that with anhydrotetracycline (aTc). Even though we later used degradation tag to accelerate the degradation process of Cre and to decrease the expression level of Cre, considering the fact that the tet operon is less prone to leakage and that using merely lac operon to control the expression of RT and Cre may cause unexpected problems, we still used different operons to control the expression of RT and Cre. This setup will be considered in the model in Part III.
There is still something that is not well explained in our current model. The final percentage of recombined Ptarget is around 1.5%. The unit of the substance is molecules, so it means there is 0.075 recombined Ptarget in one cell, which is unrealistic. This problem reflects that converting the unit of substance into molecule when doing deterministic modelling cannot offer a precise description of the system’s status.

We then used stochastic modelling techniques to determine whether and how many recombined Ptargets will show up in one replication cycle of E. coli.
Figure 4. Cre recombination (deterministic).

a) Schematic diagram of the model. Ps, Un-recombined Ptarget. Pp, Recombined Ptarget. b-c) Recombination when Cre is expressed under Lac operon. Dynamics of the percentage of un-recombined/ recombined Ptarget among all Ptargets is shown in b. Horizontal axis shows the length of time (8 hours, corresponding to R-Evolution’s function period). The distribution of the percentage of substances at the steady-state is shown in c. d) Yield of recombined Ptarget at different initial number of cDNA and Cre. The yield of recombined Ptarget is calculated as the percentage of recombined Ptarget among all Ptargets. The horizontal white line corresponds to current situation where the initial number of cDNA is 5 molecules in one E.coli. e) Yield of recombined Ptarget at different initial number of Cre when initial number of cDNA is 5 molecules. f-g) Recombination when Cre is expressed under different operon. Dynamics of the percentage of un-recombined/ recombined Ptarget among all Ptargets is shown in f. The distribution of the percentage of substances at the steady-state is shown in g.

Part II: Stochastic model to compute times of occurrence of recombined Ptarget
We use Gillespie algorithm in stochastic modelling. The procedure of this algorithm is shown as follows in the form of pseudocode:

Step 1: Initialize the reaction system at \(t=0\) with rate constants \(c_1, c_2, ......, c_v\) as initial numbers of molecules \(x_1, x_2, ......, x_u\) corresponding to \(v\) reactions and \(u\) sustances (both reactants and products) involved in the reaction system.

Step 2: For each \(i=1,2,......,v\), calculate the hazard for the \(i^{th}\) type of reaction, denoted as \(h_i(x,c_i)\) based on current substance number \(x\).

Step 3: Calculate the combined reaction hazard \(h_0(x,c)=\Sigma_{i=1}^{v}h_i(x,c_i)\).

Step 4: Simulate the time to the next reaction, \(t^{'}\) , which is a random quantity that obeys exponential distribution with parameter \(\lambda\).

Step 5: Put \(t:=t+t^{'}\).

Step 6: Simulate the reaction index, \(j\). The probability that the \(j^{th}\) reaction can occur is \(\frac{h_i(x,c_i)}{h_0(x,c)}, i=1,2,......,v.\).

Step 7: Update \(x\) according to reaction \(j\), which means putting \(x:=x+S^{(j)}\), where \(S^{(j)}\) denotes the \(j^{th}\) colomn of the stoichiometry matrix \(S\). The \(j^{th}\) column of denotes the change in substance number caused by the \(j^{th}\) reaction.

Step 8: Record time \(t\) and current substance number \(x\).

Step 9: If \(t\ \text{<}\ T_{max}\), return to step 2. \(T_{max}\) corresponds to the maximum duration of the reaction set by the user.

Step 10: Plot the result to see the dynamic of the quantity of the substance that you are interested in.

Although the algorithm is rather simple, basic mathematical skills is required to understand its theoretical basis. You may consult the book written by Wilkinson and Darren J. for a thorough understanding. The result is shown in Fig. 5.

The result demonstrates that recombined Ptargets do occur and two rounds of reverse transcription and recombination can take place in one replication cycle of E. coli (1200 s) (Fig 5a). On the contrary, no recombined will come out within that period if the initial cDNA is 5 molecules and initial Cre is 3228 molecules (Fig 5b). This again demonstrates the necessity of putting RT and Cre under different induction setups. The fluctuation of the number of recombined Ptarget is due to the backward reaction that Cre can rebind with recombined and reverting the action, making it not counted as recombined Ptarget by the algorithm.
Figure 5. Cre recombination (stochastic).

Horizontal axis shows the length of time (20min, corresponding to the duration of 1 E.coli replication cycle). Vertical axis demonstrates the number of recombined Ptarget. The initial number of Cre is 2 molecules in a, 3228 molecules in b.

Part III: Deterministic model to determine optimal induction time
To ensure the evolved protein is encoded by the mutated GOI sequence that is recombined into Ptarget, we decided to use degradation tag to accelerate the degradation process of Cre. This design would make Cre only function when inducer is in the system, thus allowing stringent control of the protein. However, we then face the problem of how to select the optimal degradation tag. Empirically, to minimize the duration of recombination, we tend to choose degradation tags with higher efficiency, but extremely high degradation rate will also reduce the yield of recombined Ptarget, leading to decreased library size. Also, it is impractical for researchers to do experiments to test these degradation tags one by one. For these reasons, we are going to use models to find out the optimal degradation tag that should be added to Cre based on the average yield of recombined Ptarget at the end of R-Evolution functioning period (8 hours).

We intend to use the models described in Part I, combined with aTc induction model proposed by Steel et al. to compute the yield of recombined Ptarget under different degradation rate of Cre (the reason why Tet operon is used has been elaborated in Part I; the schematic diagram of this process is shown in Fig 6a. The ODEs describing these processes are shown as follows. Details of the substance names, parameter names and chemical equations we used can be found in the appendix.
Although the setup in Part I successfully provided us with a clear insight into the reactions and dynamic changes of substances that underlie our mutagenesis system, the simplification that the steady-state substance concentrations of previous models can be used as inputs for subsequent models doesn’t match real reaction situation. For example, when Cre is expressed, it can immediately bind with cDNA and initiate recombination. This fact contradicts with our model assumption that recombination only takes place after both Cre and cDNA has reached their steady-state concentration.

To overcome this issue, we decided to combine all three minor models together and calculate the expected output.

As a result of the impreciseness of the basic assumption of the models in part I, we only gave a qualitative conclusion that the amount of RT and Cre should be different. Here we need to quantify how Cre degradation rate and steady-state concentration affects the yield of recombined Ptarget. That’s why we employed deterministic model here to combine the separate steps together into one and better simulate real intracellular circumstances.

By combining the models that have been talked above, we revealed the reason why the degradation tag with a moderate degradation rate, which can’t be too high or too low, should be selected (Fig 6b): under appropriate inducer concentration (20~22uM), when the degradation rate is relatively low (below 0.1 min-1), the yield of recombined Ptarget will increase according to the increase of Cre degradation rate, but when that rate is sufficiently high (above 0.1 min-1), the increase of Cre degradation rate will do harm to the yield of recombined Ptarget.

The average degradation rate acquired from literature is 0.2 min-1(Nikos et al.) and the degradation rate of Cre when tagged with the most efficient degradation tag is 0.69 min-1. Within this range of degradation rate, the maximum yield of recombined Ptarget will decrease according to the increase of Cre degradation efficiency (Fig 6c). So we decided to use the least efficient degradation tag.

We also revealed the dynamic change of the recombined Ptarget. It will continuously accumulate within Cre function period (Fig 6d). However, the concentration remains to be low within that period, due to Cre degradation (Fig 6e).

Finally, there is another interesting phenomenon that is worth mentioning. From Fig 6b and Fig 6c, we can find that for each degradation tag rate greater than 0.2 min-1, there exits a range of aTc dosage that can make the yield of recombined relatively big. Also, decreased degradation efficiency enlarges that range. This discovery provides us with another reason for using less efficient degradation tag in that it can increase the robustness of our mutagenesis system by decreasing its sensitivity to the change of inducer dosage.
Figure 6. Whole process simulation considering Cre degradation tag.

a) Chemical reactions involved in Cre induced expression. b) Yield of recombined Ptarget at different Cre degradation rate and aTc dosage. The white line on the left corresponding to the case where the degradation rate of Cre is 0.2 min-1. The white line on the left corresponding to the case where the degradation rate of Cre is 0.69 min-1. c) Yield of recombined Ptarget at different Cre degradation rate and aTc dosage (3D plot). The range of Cre degradation rate is 0.2~0.69 min-1. d) Dynamics of yield of recombined Ptarget at the Cre degradation rate of 0.2 min-1 and the initial 22uM aTc dosage. e) Dynamics of yield of recombined Ptarget at the Cre degradation rate of 0.2 min-1 and the initial 22uM aTc dosage.

In the deterministic model, we combined the three minor models proposed previously and assessed the mutagenesis system in whole. Through this addition, we achieved a better simulation of the real intracellular reactions and answered the question of when Cre should be induced for the highest level of recombination efficiency to be obtained.
Complementary materials
You can access the Matlab codes of our models from our Github repository.

Please consult the appendix file for a clearer understanding of the formulation of the model.

References