Skip to content

Instantly share code, notes, and snippets.

@celskeggs
Last active December 7, 2019 01:04
Show Gist options
  • Save celskeggs/99ba814e9b0e1896081d3c02f3423978 to your computer and use it in GitHub Desktop.
Save celskeggs/99ba814e9b0e1896081d3c02f3423978 to your computer and use it in GitHub Desktop.
1-bit Serial Register supplementary data

Supplementary material for: https://www.youtube.com/watch?v=yfWoBryKI9U

Run generator with:

$ python3 generate_matlab.py

Download dependencies:

$ wget https://web.mit.edu/20.305/www/v5/Rate.m
$ wget https://web.mit.edu/20.305/www/v5/Pulse.m
$ wget https://web.mit.edu/20.305/www/v5/Part.m
$ wget https://web.mit.edu/20.305/www/v5/Const.m
$ wget https://web.mit.edu/20.305/www/v5/Compositor.m
$ wget https://web.mit.edu/20.305/www/v5/BioSystem.m
$ sha256sum *.m
9819ed5940013e41e7f0b7375d98d4e73b1a5d4eafb231a2ba3405bf252f33ce  BioSystem.m
01170babb1f1946e23f123e1c798c108017f2da3429cf685eee68842aa5f75cb  Compositor.m
10f4edfcea41ad0fb21f6eafdeaf6aae33ab723d021236f68004ee6f4db7d766  Const.m
bc1624e39f3c613e81f2677f364c26c0ae1196d1e37a0829561e25666466a5df  Part.m
73e31cfab8e07004620ec81e0add51b6e41733bb27e8dadc0af5d1a87885fddc  Pulse.m
db9fb4579b1322795c89ab04ba9e2fd0d11313521ecb223a981c1153b7cc2240  Rate.m

Run project in MATLAB with:

>> project
#!/usr/bin/env python3
NONRATE = "*-^()/+"
def enum_vars_simple(x):
all_vars = set()
for nr in NONRATE:
x = x.replace(nr, " ")
names = x.split()
for name in names:
if name.replace(".", "").isdigit():
continue
assert name.replace("_", "").isalnum(), "invalid name: %s" % name
all_vars.add(name)
return all_vars
def enum_vars(rates):
all_vars = set(rates.keys())
for x in rates.values():
all_vars.update(enum_vars_simple(x))
return all_vars
class System:
def __init__(self, simlength, ylim=None):
self.constants = {}
self.initials = {}
self.vars = set()
self.parts = {}
self.pulses = []
self.plotcolors = {}
self.plots = []
self.simlength = simlength
self.ylim = ylim
def add_const(self, name, value):
assert name not in self.constants
assert name not in self.vars
assert name not in self.initials
self.constants[name] = value
def add_part(self, desc, rates):
assert set(rates.keys()).isdisjoint(set(self.constants.keys()))
all_vars = enum_vars(rates).difference(set(self.constants.keys()))
self.vars = self.vars.union(all_vars)
assert desc not in self.parts, "already found: %s" % desc
self.parts[desc] = {var: str(rates.get(var, "0")) for var in all_vars}
def express(self, rate, *names):
for name in names:
assert name and name not in self.constants
self.add_part("expression of %s" % name, {name: rate})
def activated(self, *names, basal=None, maximal=None, RX=()):
assert basal is not None and maximal is not None
for name in sorted(set(names)):
mul = names.count(name)
hillTerms = []
for RI in RX:
if type(RI) == tuple:
R, I = RI
hillI = "{I}^hill_{I} / (threshold_{I}^hill_{I} + {I}^hill_{I})".format(I=I)
else:
assert type(RI) == str
R = RI
hillI = "1"
hillTerms += [" * (threshold_{R}^hill_{R} / (threshold_{R}^hill_{R} + ({R} * ({HI}))^hill_{R}))".format(R=R, HI=hillI)]
self.add_part("expression of %s" % name, {name: "{mul} * ({basal} + ({maximal} - {basal}){hill})".format(mul=mul, basal=basal, maximal=maximal, hill="".join(hillTerms))})
def dimer(self, bind_rate, unbind_rate, half_a, half_b, whole):
dimerization = "{bind} * {half_a} * {half_b}".format(bind=bind_rate, half_a=half_a, half_b=half_b)
dedimerization = "{unbind} * {whole}".format(unbind=unbind_rate, whole=whole)
if half_a == half_b:
self.add_part("homodimerization of %s" % (half_a), {
half_a: "-2*(%s)" % dimerization,
whole: dimerization,
})
self.add_part("dedimerization of %s" % whole, {
half_a: "2*(%s)" % dedimerization,
whole: "-(%s)" % dedimerization,
})
else:
self.add_part("dimerization of %s and %s" % (half_a, half_b), {
half_a: "-(%s)" % dimerization,
half_b: "-(%s)" % dimerization,
whole: dimerization,
})
self.add_part("dedimerization of %s" % whole, {
half_a: dedimerization,
half_b: dedimerization,
whole: "-(%s)" % dedimerization,
})
def degrade(self, rate, *names):
assert rate in self.constants
for name in names:
assert name not in self.constants
self.add_part("degradation of %s" % name, {name: "- %s * %s" % (rate, name)})
def plot(self, *lines, colors=()):
self.plotcolors.update({line: color for line, color in zip(lines, colors)})
self.plots.append(sorted(lines))
def generate(self, wl):
wl("% generated automatically by cela's gen.py")
wl("clear all;")
wl("sys = BioSystem();")
wl("")
wl("% constants")
for k, v in sorted(self.constants.items()):
wl("sys.AddConstant(Const('%s', %f));" % (k, v))
wl("")
wl("% compositors")
for k in sorted(set(self.initials.keys()).union(self.vars)):
wl("d%sdt = sys.AddCompositor('%s', %f);" % (k, k, self.initials.get(k, 0)))
wl("")
wl("% parts")
for desc, rates in self.parts.items():
order = sorted(rates.keys())
compositors = ["d%sdt" % ref for ref in order]
wl("sys.AddPart(Part('%s', [%s], [" % (desc, " ".join(compositors)))
for ref in order:
wl(" Rate('%s')" % rates[ref])
wl(" ]));")
wl("")
wl("% simulate")
if self.pulses:
wl("[T, Y] = sys.run_pulses([")
pulses = sorted(self.pulses, key=lambda p: p[0])
if pulses[0][0] > 0:
somename = pulses[0][1]
pulses.insert(0, (0, somename, self.initials.get(somename, 0)))
for time, var, concentration in pulses:
assert var in self.vars or var in self.initials
wl(" Pulse(%f, '%s', %f)" % (time, var, concentration))
wl(" Pulse(%f, '', 0) %% stop simulation" % (self.simlength))
wl("]);")
else:
wl("[T, Y] = sys.run([0 %f]);" % (self.simlength))
if self.plots:
wl("")
wl("% plot")
wl("figure(1);");
for i, plot in enumerate(self.plots, 1):
wl("subplot(%d,1,%d);" % (len(self.plots), i))
wl("plot(...");
for line in plot:
wl(" T, Y(:, sys.CompositorIndex('%s')), '%s', ..." % (line, self.plotcolors.get(line, "")))
wl(" 'LineWidth', 2);")
wl("legend(%s, 'Location', 'northeast')" % ", ".join("'%s'" % line for line in plot))
wl("xlabel('Time (min)');")
wl("ylabel('Concentration (nM)');")
if self.ylim is not None:
wl("ylim([0 %f]);" % self.ylim)
def write_to(self, filename):
with open(filename, "w") as f:
self.generate(lambda x: f.write(x + "\n"))
def writer(f):
return lambda filename: f().write_to(filename)
tetRs = { # threshold, hill, max, basal
'TetR': (0.1, 2.7, 24, 0.2),
'QacR': (0.5, 1.4, 21, 0.2),
'IcaR': (0.4, 1.8, 13, 0.4),
'AmeR': (0.5, 1.4, 17, 1.7),
'ScbR': (0.2, 2.6, 5, 0.6),
'LmrA': (1.2, 3.1, 70, 1.1),
'AmtR': (0.2, 1.8, 9, 0.3),
'SmcR': (0.1, 2, 13, 2.1),
'McbR': (0.4, 1.6, 16, 1.1),
'BetI': (0.2, 2.4, 13, 0.4),
'SrpR': (0.3, 3.2, 25, 0.1),
'Orf2': (0.4, 6.1, 14, 0.2),
'BM3R1': (0.6, 4.5, 3, 0.1),
'TarA': (0.1, 1.8, 13, 0.2),
'PhlF': (0.4, 4.5, 16, 0.1),
'ButR': (1.3, 2.4, 12, 1.8),
'PsrA': (0.4, 2, 20, 0.5),
'HapR': (0.2, 1.4, 10, 0.9),
'HlyIIR': (0.5, 2.7, 17, 0.3),
'LitR': (0.1, 1.9, 16, 0.5),
}
@writer
def project():
sys = System(simlength=8000)
REU = 1.2
THR = 1.0 / 0.0165
patched = True
toggle1a = "TetR"
toggle1b = "BM3R1"
toggle2a = "SrpR"
toggle2b = "PhlF"
middle = "QacR"
needed = [toggle1a, toggle1b, toggle2a, toggle2b, middle]
for name in needed:
n_threshold, n_hill, n_max, n_basal = tetRs[name]
sys.constants["basal_" + name] = n_basal * REU
sys.constants["max_" + name] = n_max * REU
sys.constants["threshold_" + name] = n_threshold * REU * THR
sys.constants["hill_" + name] = n_hill
sys.constants["threshold_IPTG"] = 0.2926
sys.constants["threshold_INDNX"] = 0.2926
sys.constants["hill_IPTG"] = 2.0
sys.constants["hill_INDNX"] = 2.0
sys.constants["basal_ISV94A"] = 0.8 * REU
sys.constants["max_ISV94A"] = 5.5 * REU
sys.constants["threshold_ISV94A"] = 0.2
sys.constants["threshold_ISNX"] = 0.2
sys.constants["hill_ISV94A"] = 2.0
sys.constants["hill_ISNX"] = 2.0
sys.constants["hill_R434hR43422"] = 1.0
sys.constants["hill_R434hR434NX"] = 1.0
sys.constants["threshold_R434hR434NX"] = 15 if patched else 5.0
sys.constants["threshold_R434hR43422"] = 25 if patched else 5.0
sys.constants["bind_R434"] = 0.2
sys.constants["unbind_R434"] = 2
sys.constants["constitutive"] = 1.0 * REU
sys.constants["degrade_protein"] = 0.0165
sys.constants["degrade_pdt3"] = 0.0825
sys.constants["degrade_small_molecule"] = 0.0165
IPTG, INDNX = "IPTG", "INDNX"
ISV94A, ISNX = "ISV94A", "ISNX"
BM3R1, ICAR, LMRA, PHLF, SRPR = toggle1a, toggle1b, middle, toggle2a, toggle2b
R434, R43422, R434NX = "R434", "R43422", "R434NX"
H_R4_R4 = R434 + "h" + R434
H_R4_22 = R434 + "h" + R43422
H_R4_NX = R434 + "h" + R434NX
H_22_22 = R43422 + "h" + R43422
H_22_NX = R43422 + "h" + R434NX
H_NX_NX = R434NX + "h" + R434NX
YFP, RFP, GFP = "YFP", "RFP", "GFP"
sys.degrade("degrade_protein", YFP, RFP, GFP, ISV94A, ISNX, R434, R43422, R434NX, H_R4_R4, H_R4_22, H_R4_NX, H_22_22, H_22_NX, H_NX_NX, SRPR, PHLF, BM3R1, ICAR, LMRA)
sys.degrade("degrade_small_molecule", IPTG, INDNX)
# pulse concentration
pc = 10 # nM
sys.pulses.append((500, IPTG, pc))
sys.pulses.append((2000, IPTG, pc))
sys.pulses.append((3500, INDNX, pc))
sys.pulses.append((5000, INDNX, pc))
sys.pulses.append((6500, IPTG, pc))
sys.express("constitutive", ISV94A, ISNX)
sys.activated(ICAR, R43422, YFP, basal="basal_" + BM3R1, maximal="max_" + BM3R1, RX=[(ISNX, INDNX), BM3R1])
sys.activated(BM3R1, R434NX, basal="basal_" + ICAR, maximal="max_" + ICAR, RX=[(ISV94A, IPTG), ICAR])
sys.activated(LMRA, R434, basal="basal_" + ISV94A, maximal="max_" + ISV94A, RX=[(ISNX, INDNX), (ISV94A, IPTG)])
sys.dimer("bind_R434", "unbind_R434", R434, R434, H_R4_R4)
sys.dimer("bind_R434", "unbind_R434", R434, R43422, H_R4_22)
sys.dimer("bind_R434", "unbind_R434", R434, R434NX, H_R4_NX)
sys.dimer("bind_R434", "unbind_R434", R43422, R43422, H_22_22)
sys.dimer("bind_R434", "unbind_R434", R43422, R434NX, H_22_NX)
sys.dimer("bind_R434", "unbind_R434", R434NX, R434NX, H_NX_NX)
if patched:
sys.activated(SRPR, basal="basal_" + PHLF, maximal="max_" + PHLF, RX=[H_R4_NX, H_R4_NX, PHLF])
sys.activated(PHLF, basal="basal_" + SRPR, maximal="max_" + SRPR, RX=[H_R4_22, H_R4_22, SRPR])
else:
sys.activated(SRPR, basal="basal_" + PHLF, maximal="max_" + PHLF, RX=[H_R4_NX, PHLF])
sys.activated(PHLF, basal="basal_" + SRPR, maximal="max_" + SRPR, RX=[H_R4_22, SRPR])
sys.activated(GFP, basal="basal_" + PHLF, maximal="max_" + PHLF, RX=[LMRA, PHLF])
if patched:
sys.activated(RFP, RFP, basal="basal_" + SRPR, maximal="max_" + SRPR, RX=[LMRA, SRPR])
else:
sys.activated(RFP, basal="basal_" + SRPR, maximal="max_" + SRPR, RX=[LMRA, SRPR])
sys.plot(IPTG, INDNX, colors=('g', 'r'))
sys.plot(LMRA, ICAR, BM3R1, colors=('m', 'g', 'r'))
sys.plot(R434, H_R4_22, H_R4_NX, colors=('m', 'g', 'r'))
sys.plot(PHLF, SRPR, colors=('r', 'g'))
sys.plot(YFP, GFP, RFP, colors=('y', 'g', 'r'))
expected_vars = {ISV94A, ISNX, LMRA, BM3R1, ICAR, PHLF, SRPR, R434, R43422, R434NX, YFP, RFP, GFP, IPTG, INDNX, H_R4_R4, H_R4_22, H_R4_NX, H_22_22, H_22_NX, H_NX_NX}
print("all variables:", sys.vars)
unexpected = sys.vars.difference(expected_vars)
if unexpected:
raise Exception("bad variables: %s" % " ".join(unexpected))
return sys
if __name__ == "__main__":
project("project.m")
print("generated")
% generated automatically by cela's gen.py
clear all;
sys = BioSystem();
% constants
sys.AddConstant(Const('basal_BM3R1', 0.120000));
sys.AddConstant(Const('basal_ISV94A', 0.960000));
sys.AddConstant(Const('basal_PhlF', 0.120000));
sys.AddConstant(Const('basal_QacR', 0.240000));
sys.AddConstant(Const('basal_SrpR', 0.120000));
sys.AddConstant(Const('basal_TetR', 0.240000));
sys.AddConstant(Const('bind_R434', 0.200000));
sys.AddConstant(Const('constitutive', 1.200000));
sys.AddConstant(Const('degrade_pdt3', 0.082500));
sys.AddConstant(Const('degrade_protein', 0.016500));
sys.AddConstant(Const('degrade_small_molecule', 0.016500));
sys.AddConstant(Const('hill_BM3R1', 4.500000));
sys.AddConstant(Const('hill_INDNX', 2.000000));
sys.AddConstant(Const('hill_IPTG', 2.000000));
sys.AddConstant(Const('hill_ISNX', 2.000000));
sys.AddConstant(Const('hill_ISV94A', 2.000000));
sys.AddConstant(Const('hill_PhlF', 4.500000));
sys.AddConstant(Const('hill_QacR', 1.400000));
sys.AddConstant(Const('hill_R434hR43422', 1.000000));
sys.AddConstant(Const('hill_R434hR434NX', 1.000000));
sys.AddConstant(Const('hill_SrpR', 3.200000));
sys.AddConstant(Const('hill_TetR', 2.700000));
sys.AddConstant(Const('max_BM3R1', 3.600000));
sys.AddConstant(Const('max_ISV94A', 6.600000));
sys.AddConstant(Const('max_PhlF', 19.200000));
sys.AddConstant(Const('max_QacR', 25.200000));
sys.AddConstant(Const('max_SrpR', 30.000000));
sys.AddConstant(Const('max_TetR', 28.800000));
sys.AddConstant(Const('threshold_BM3R1', 43.636364));
sys.AddConstant(Const('threshold_INDNX', 0.292600));
sys.AddConstant(Const('threshold_IPTG', 0.292600));
sys.AddConstant(Const('threshold_ISNX', 0.200000));
sys.AddConstant(Const('threshold_ISV94A', 0.200000));
sys.AddConstant(Const('threshold_PhlF', 29.090909));
sys.AddConstant(Const('threshold_QacR', 36.363636));
sys.AddConstant(Const('threshold_R434hR43422', 25.000000));
sys.AddConstant(Const('threshold_R434hR434NX', 15.000000));
sys.AddConstant(Const('threshold_SrpR', 21.818182));
sys.AddConstant(Const('threshold_TetR', 7.272727));
sys.AddConstant(Const('unbind_R434', 2.000000));
% compositors
dBM3R1dt = sys.AddCompositor('BM3R1', 0.000000);
dGFPdt = sys.AddCompositor('GFP', 0.000000);
dINDNXdt = sys.AddCompositor('INDNX', 0.000000);
dIPTGdt = sys.AddCompositor('IPTG', 0.000000);
dISNXdt = sys.AddCompositor('ISNX', 0.000000);
dISV94Adt = sys.AddCompositor('ISV94A', 0.000000);
dPhlFdt = sys.AddCompositor('PhlF', 0.000000);
dQacRdt = sys.AddCompositor('QacR', 0.000000);
dR434dt = sys.AddCompositor('R434', 0.000000);
dR43422dt = sys.AddCompositor('R43422', 0.000000);
dR43422hR43422dt = sys.AddCompositor('R43422hR43422', 0.000000);
dR43422hR434NXdt = sys.AddCompositor('R43422hR434NX', 0.000000);
dR434NXdt = sys.AddCompositor('R434NX', 0.000000);
dR434NXhR434NXdt = sys.AddCompositor('R434NXhR434NX', 0.000000);
dR434hR434dt = sys.AddCompositor('R434hR434', 0.000000);
dR434hR43422dt = sys.AddCompositor('R434hR43422', 0.000000);
dR434hR434NXdt = sys.AddCompositor('R434hR434NX', 0.000000);
dRFPdt = sys.AddCompositor('RFP', 0.000000);
dSrpRdt = sys.AddCompositor('SrpR', 0.000000);
dTetRdt = sys.AddCompositor('TetR', 0.000000);
dYFPdt = sys.AddCompositor('YFP', 0.000000);
% parts
sys.AddPart(Part('degradation of YFP', [dYFPdt], [
Rate('- degrade_protein * YFP')
]));
sys.AddPart(Part('degradation of RFP', [dRFPdt], [
Rate('- degrade_protein * RFP')
]));
sys.AddPart(Part('degradation of GFP', [dGFPdt], [
Rate('- degrade_protein * GFP')
]));
sys.AddPart(Part('degradation of ISV94A', [dISV94Adt], [
Rate('- degrade_protein * ISV94A')
]));
sys.AddPart(Part('degradation of ISNX', [dISNXdt], [
Rate('- degrade_protein * ISNX')
]));
sys.AddPart(Part('degradation of R434', [dR434dt], [
Rate('- degrade_protein * R434')
]));
sys.AddPart(Part('degradation of R43422', [dR43422dt], [
Rate('- degrade_protein * R43422')
]));
sys.AddPart(Part('degradation of R434NX', [dR434NXdt], [
Rate('- degrade_protein * R434NX')
]));
sys.AddPart(Part('degradation of R434hR434', [dR434hR434dt], [
Rate('- degrade_protein * R434hR434')
]));
sys.AddPart(Part('degradation of R434hR43422', [dR434hR43422dt], [
Rate('- degrade_protein * R434hR43422')
]));
sys.AddPart(Part('degradation of R434hR434NX', [dR434hR434NXdt], [
Rate('- degrade_protein * R434hR434NX')
]));
sys.AddPart(Part('degradation of R43422hR43422', [dR43422hR43422dt], [
Rate('- degrade_protein * R43422hR43422')
]));
sys.AddPart(Part('degradation of R43422hR434NX', [dR43422hR434NXdt], [
Rate('- degrade_protein * R43422hR434NX')
]));
sys.AddPart(Part('degradation of R434NXhR434NX', [dR434NXhR434NXdt], [
Rate('- degrade_protein * R434NXhR434NX')
]));
sys.AddPart(Part('degradation of PhlF', [dPhlFdt], [
Rate('- degrade_protein * PhlF')
]));
sys.AddPart(Part('degradation of SrpR', [dSrpRdt], [
Rate('- degrade_protein * SrpR')
]));
sys.AddPart(Part('degradation of TetR', [dTetRdt], [
Rate('- degrade_protein * TetR')
]));
sys.AddPart(Part('degradation of BM3R1', [dBM3R1dt], [
Rate('- degrade_protein * BM3R1')
]));
sys.AddPart(Part('degradation of QacR', [dQacRdt], [
Rate('- degrade_protein * QacR')
]));
sys.AddPart(Part('degradation of IPTG', [dIPTGdt], [
Rate('- degrade_small_molecule * IPTG')
]));
sys.AddPart(Part('degradation of INDNX', [dINDNXdt], [
Rate('- degrade_small_molecule * INDNX')
]));
sys.AddPart(Part('expression of ISV94A', [dISV94Adt], [
Rate('constitutive')
]));
sys.AddPart(Part('expression of ISNX', [dISNXdt], [
Rate('constitutive')
]));
sys.AddPart(Part('expression of BM3R1', [dBM3R1dt dINDNXdt dISNXdt dTetRdt], [
Rate('1 * (basal_TetR + (max_TetR - basal_TetR) * (threshold_ISNX^hill_ISNX / (threshold_ISNX^hill_ISNX + (ISNX * (INDNX^hill_INDNX / (threshold_INDNX^hill_INDNX + INDNX^hill_INDNX)))^hill_ISNX)) * (threshold_TetR^hill_TetR / (threshold_TetR^hill_TetR + (TetR * (1))^hill_TetR)))')
Rate('0')
Rate('0')
Rate('0')
]));
sys.AddPart(Part('expression of R43422', [dINDNXdt dISNXdt dR43422dt dTetRdt], [
Rate('0')
Rate('0')
Rate('1 * (basal_TetR + (max_TetR - basal_TetR) * (threshold_ISNX^hill_ISNX / (threshold_ISNX^hill_ISNX + (ISNX * (INDNX^hill_INDNX / (threshold_INDNX^hill_INDNX + INDNX^hill_INDNX)))^hill_ISNX)) * (threshold_TetR^hill_TetR / (threshold_TetR^hill_TetR + (TetR * (1))^hill_TetR)))')
Rate('0')
]));
sys.AddPart(Part('expression of YFP', [dINDNXdt dISNXdt dTetRdt dYFPdt], [
Rate('0')
Rate('0')
Rate('0')
Rate('1 * (basal_TetR + (max_TetR - basal_TetR) * (threshold_ISNX^hill_ISNX / (threshold_ISNX^hill_ISNX + (ISNX * (INDNX^hill_INDNX / (threshold_INDNX^hill_INDNX + INDNX^hill_INDNX)))^hill_ISNX)) * (threshold_TetR^hill_TetR / (threshold_TetR^hill_TetR + (TetR * (1))^hill_TetR)))')
]));
sys.AddPart(Part('expression of R434NX', [dBM3R1dt dIPTGdt dISV94Adt dR434NXdt], [
Rate('0')
Rate('0')
Rate('0')
Rate('1 * (basal_BM3R1 + (max_BM3R1 - basal_BM3R1) * (threshold_ISV94A^hill_ISV94A / (threshold_ISV94A^hill_ISV94A + (ISV94A * (IPTG^hill_IPTG / (threshold_IPTG^hill_IPTG + IPTG^hill_IPTG)))^hill_ISV94A)) * (threshold_BM3R1^hill_BM3R1 / (threshold_BM3R1^hill_BM3R1 + (BM3R1 * (1))^hill_BM3R1)))')
]));
sys.AddPart(Part('expression of TetR', [dBM3R1dt dIPTGdt dISV94Adt dTetRdt], [
Rate('0')
Rate('0')
Rate('0')
Rate('1 * (basal_BM3R1 + (max_BM3R1 - basal_BM3R1) * (threshold_ISV94A^hill_ISV94A / (threshold_ISV94A^hill_ISV94A + (ISV94A * (IPTG^hill_IPTG / (threshold_IPTG^hill_IPTG + IPTG^hill_IPTG)))^hill_ISV94A)) * (threshold_BM3R1^hill_BM3R1 / (threshold_BM3R1^hill_BM3R1 + (BM3R1 * (1))^hill_BM3R1)))')
]));
sys.AddPart(Part('expression of QacR', [dINDNXdt dIPTGdt dISNXdt dISV94Adt dQacRdt], [
Rate('0')
Rate('0')
Rate('0')
Rate('0')
Rate('1 * (basal_ISV94A + (max_ISV94A - basal_ISV94A) * (threshold_ISNX^hill_ISNX / (threshold_ISNX^hill_ISNX + (ISNX * (INDNX^hill_INDNX / (threshold_INDNX^hill_INDNX + INDNX^hill_INDNX)))^hill_ISNX)) * (threshold_ISV94A^hill_ISV94A / (threshold_ISV94A^hill_ISV94A + (ISV94A * (IPTG^hill_IPTG / (threshold_IPTG^hill_IPTG + IPTG^hill_IPTG)))^hill_ISV94A)))')
]));
sys.AddPart(Part('expression of R434', [dINDNXdt dIPTGdt dISNXdt dISV94Adt dR434dt], [
Rate('0')
Rate('0')
Rate('0')
Rate('0')
Rate('1 * (basal_ISV94A + (max_ISV94A - basal_ISV94A) * (threshold_ISNX^hill_ISNX / (threshold_ISNX^hill_ISNX + (ISNX * (INDNX^hill_INDNX / (threshold_INDNX^hill_INDNX + INDNX^hill_INDNX)))^hill_ISNX)) * (threshold_ISV94A^hill_ISV94A / (threshold_ISV94A^hill_ISV94A + (ISV94A * (IPTG^hill_IPTG / (threshold_IPTG^hill_IPTG + IPTG^hill_IPTG)))^hill_ISV94A)))')
]));
sys.AddPart(Part('homodimerization of R434', [dR434dt dR434hR434dt], [
Rate('-2*(bind_R434 * R434 * R434)')
Rate('bind_R434 * R434 * R434')
]));
sys.AddPart(Part('dedimerization of R434hR434', [dR434dt dR434hR434dt], [
Rate('2*(unbind_R434 * R434hR434)')
Rate('-(unbind_R434 * R434hR434)')
]));
sys.AddPart(Part('dimerization of R434 and R43422', [dR434dt dR43422dt dR434hR43422dt], [
Rate('-(bind_R434 * R434 * R43422)')
Rate('-(bind_R434 * R434 * R43422)')
Rate('bind_R434 * R434 * R43422')
]));
sys.AddPart(Part('dedimerization of R434hR43422', [dR434dt dR43422dt dR434hR43422dt], [
Rate('unbind_R434 * R434hR43422')
Rate('unbind_R434 * R434hR43422')
Rate('-(unbind_R434 * R434hR43422)')
]));
sys.AddPart(Part('dimerization of R434 and R434NX', [dR434dt dR434NXdt dR434hR434NXdt], [
Rate('-(bind_R434 * R434 * R434NX)')
Rate('-(bind_R434 * R434 * R434NX)')
Rate('bind_R434 * R434 * R434NX')
]));
sys.AddPart(Part('dedimerization of R434hR434NX', [dR434dt dR434NXdt dR434hR434NXdt], [
Rate('unbind_R434 * R434hR434NX')
Rate('unbind_R434 * R434hR434NX')
Rate('-(unbind_R434 * R434hR434NX)')
]));
sys.AddPart(Part('homodimerization of R43422', [dR43422dt dR43422hR43422dt], [
Rate('-2*(bind_R434 * R43422 * R43422)')
Rate('bind_R434 * R43422 * R43422')
]));
sys.AddPart(Part('dedimerization of R43422hR43422', [dR43422dt dR43422hR43422dt], [
Rate('2*(unbind_R434 * R43422hR43422)')
Rate('-(unbind_R434 * R43422hR43422)')
]));
sys.AddPart(Part('dimerization of R43422 and R434NX', [dR43422dt dR43422hR434NXdt dR434NXdt], [
Rate('-(bind_R434 * R43422 * R434NX)')
Rate('bind_R434 * R43422 * R434NX')
Rate('-(bind_R434 * R43422 * R434NX)')
]));
sys.AddPart(Part('dedimerization of R43422hR434NX', [dR43422dt dR43422hR434NXdt dR434NXdt], [
Rate('unbind_R434 * R43422hR434NX')
Rate('-(unbind_R434 * R43422hR434NX)')
Rate('unbind_R434 * R43422hR434NX')
]));
sys.AddPart(Part('homodimerization of R434NX', [dR434NXdt dR434NXhR434NXdt], [
Rate('-2*(bind_R434 * R434NX * R434NX)')
Rate('bind_R434 * R434NX * R434NX')
]));
sys.AddPart(Part('dedimerization of R434NXhR434NX', [dR434NXdt dR434NXhR434NXdt], [
Rate('2*(unbind_R434 * R434NXhR434NX)')
Rate('-(unbind_R434 * R434NXhR434NX)')
]));
sys.AddPart(Part('expression of PhlF', [dPhlFdt dR434hR434NXdt dSrpRdt], [
Rate('1 * (basal_SrpR + (max_SrpR - basal_SrpR) * (threshold_R434hR434NX^hill_R434hR434NX / (threshold_R434hR434NX^hill_R434hR434NX + (R434hR434NX * (1))^hill_R434hR434NX)) * (threshold_R434hR434NX^hill_R434hR434NX / (threshold_R434hR434NX^hill_R434hR434NX + (R434hR434NX * (1))^hill_R434hR434NX)) * (threshold_SrpR^hill_SrpR / (threshold_SrpR^hill_SrpR + (SrpR * (1))^hill_SrpR)))')
Rate('0')
Rate('0')
]));
sys.AddPart(Part('expression of SrpR', [dPhlFdt dR434hR43422dt dSrpRdt], [
Rate('0')
Rate('0')
Rate('1 * (basal_PhlF + (max_PhlF - basal_PhlF) * (threshold_R434hR43422^hill_R434hR43422 / (threshold_R434hR43422^hill_R434hR43422 + (R434hR43422 * (1))^hill_R434hR43422)) * (threshold_R434hR43422^hill_R434hR43422 / (threshold_R434hR43422^hill_R434hR43422 + (R434hR43422 * (1))^hill_R434hR43422)) * (threshold_PhlF^hill_PhlF / (threshold_PhlF^hill_PhlF + (PhlF * (1))^hill_PhlF)))')
]));
sys.AddPart(Part('expression of GFP', [dGFPdt dQacRdt dSrpRdt], [
Rate('1 * (basal_SrpR + (max_SrpR - basal_SrpR) * (threshold_QacR^hill_QacR / (threshold_QacR^hill_QacR + (QacR * (1))^hill_QacR)) * (threshold_SrpR^hill_SrpR / (threshold_SrpR^hill_SrpR + (SrpR * (1))^hill_SrpR)))')
Rate('0')
Rate('0')
]));
sys.AddPart(Part('expression of RFP', [dPhlFdt dQacRdt dRFPdt], [
Rate('0')
Rate('0')
Rate('2 * (basal_PhlF + (max_PhlF - basal_PhlF) * (threshold_QacR^hill_QacR / (threshold_QacR^hill_QacR + (QacR * (1))^hill_QacR)) * (threshold_PhlF^hill_PhlF / (threshold_PhlF^hill_PhlF + (PhlF * (1))^hill_PhlF)))')
]));
% simulate
[T, Y] = sys.run_pulses([
Pulse(0.000000, 'IPTG', 0.000000)
Pulse(500.000000, 'IPTG', 10.000000)
Pulse(2000.000000, 'IPTG', 10.000000)
Pulse(3500.000000, 'INDNX', 10.000000)
Pulse(5000.000000, 'INDNX', 10.000000)
Pulse(6500.000000, 'IPTG', 10.000000)
Pulse(8000.000000, '', 0) % stop simulation
]);
% plot
figure(1);
subplot(5,1,1);
plot(...
T, Y(:, sys.CompositorIndex('INDNX')), 'r', ...
T, Y(:, sys.CompositorIndex('IPTG')), 'g', ...
'LineWidth', 2);
legend('INDNX', 'IPTG', 'Location', 'northeast')
xlabel('Time (min)');
ylabel('Concentration (nM)');
subplot(5,1,2);
plot(...
T, Y(:, sys.CompositorIndex('BM3R1')), 'g', ...
T, Y(:, sys.CompositorIndex('QacR')), 'm', ...
T, Y(:, sys.CompositorIndex('TetR')), 'r', ...
'LineWidth', 2);
legend('BM3R1', 'QacR', 'TetR', 'Location', 'northeast')
xlabel('Time (min)');
ylabel('Concentration (nM)');
subplot(5,1,3);
plot(...
T, Y(:, sys.CompositorIndex('R434')), 'm', ...
T, Y(:, sys.CompositorIndex('R434hR43422')), 'g', ...
T, Y(:, sys.CompositorIndex('R434hR434NX')), 'r', ...
'LineWidth', 2);
legend('R434', 'R434hR43422', 'R434hR434NX', 'Location', 'northeast')
xlabel('Time (min)');
ylabel('Concentration (nM)');
subplot(5,1,4);
plot(...
T, Y(:, sys.CompositorIndex('PhlF')), 'g', ...
T, Y(:, sys.CompositorIndex('SrpR')), 'r', ...
'LineWidth', 2);
legend('PhlF', 'SrpR', 'Location', 'northeast')
xlabel('Time (min)');
ylabel('Concentration (nM)');
subplot(5,1,5);
plot(...
T, Y(:, sys.CompositorIndex('GFP')), 'g', ...
T, Y(:, sys.CompositorIndex('RFP')), 'r', ...
T, Y(:, sys.CompositorIndex('YFP')), 'y', ...
'LineWidth', 2);
legend('GFP', 'RFP', 'YFP', 'Location', 'northeast')
xlabel('Time (min)');
ylabel('Concentration (nM)');
Display the source blob
Display the rendered blob
Raw
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment