[PYTHON] Einfach zu bedienendes E-Cell 4 Intermediate
Vorbereitung
Informationen zur Installation finden Sie im separaten Artikel.
Dieser Artikel ist eine Zwischenausgabe. Es wird davon ausgegangen, dass Sie die Anfängerausgabe verstehen.
Darüber hinaus wird in diesem Artikel die Verwendung von IPython Notebook vorausgesetzt.
Wie benutzt man
Raum in der E-Zelle 4
In der Anfängerausgabe habe ich die drei Elemente Modell, Welt und Simulator in E-Cell 4 erklärt. Außerdem habe ich eine einfache Berechnung mit ode, einem normalen Differentialgleichungslöser, und gilespie, einer probabilistischen Methode, versucht.
Ich habe eine Welt erstellt, indem ich bei Berechnungen mit Ode und Gilespie eine Größe angegeben habe. Wie wird jedoch mit dem Speicherplatz in E-Cell 4 umgegangen?
from ecell4 import *
w1 = ode.ODEWorld(Real3(1, 1, 1))
w2 = gillespie.GillespieWorld(Real3(1, 1, 1))
In Ode und Gillespie haben wir in einem Würfel mit einer Seite wie im obigen Beispiel berechnet, aber in Wirklichkeit ist nur dieses Volumen von Bedeutung.
w3 = ode.ODEWorld(Real3(2, 0.5, 1)) # is almost equivalent to 'w1'
w4 = gillespie.GillespieWorld(Real3(2, 2, 0.25)) # is almost equivalent to 'w2'
Da jedoch immerhin 1 als Volumen angegeben ist, wird das gleiche Berechnungsergebnis angezeigt.
Dies erscheint in einem System, das tatsächlich gut bewegt und räumlich einheitlich ist, wie in vitro, vernünftig. Die Zellen sind jedoch eindeutig nicht räumlich einheitlich. Um die Lokalisierung dieser Moleküle zu berücksichtigen. Erfordert weltraumbasierte biochemische Berechnungen.
In E-Cell 4 können verschiedene räumliche Ausdrücke und entsprechende Berechnungstechniken verwendet werden. Im Folgenden wird zunächst die Handhabung des Raums in E-Cell 4 anhand der räumlichen Gillespie-Methode betrachtet, die eine davon ist.
Räumliche Gillespie-Methode
In E-Cell 4 ist die räumliche Gillespie-Methode [^ 1] im "Meso" -Modul enthalten. Berechnen wir zunächst mit "run_simulation", wie wir in der Anfängerausgabe gesehen haben, ohne an irgendetwas zu denken.
[^ 1]: Richtiger ist, dass es näher an der Next Subvolume-Methode liegt. Weitere Informationen finden Sie im Dokument.
import numpy
from ecell4 import *
with reaction_rules():
A + B == C | (0.01, 0.3)
y = run_simulation(numpy.linspace(0, 10, 100), {'C': 60}, solver='meso')
\frac{d\mathrm{C}}{dt}=0.01\frac{\mathrm{A}}{V}\frac{\mathrm{B}}{V}-0.3\frac{\mathrm{C}}{V}=0\\
0.01\left(60-\mathrm{C}\right)^2=0.3\mathrm{C}\times V\\
\mathrm{C}=30
Im Vergleich zu Ode und Gillespie mag es einige Zeit gedauert haben, aber Sie erhalten fast das gleiche Ergebnis. Richtig, das Meso-Modul führt fast die gleiche Berechnung wie Gilespie ohne räumliche Einstellungen durch. Lassen Sie uns 'run_simulation' erweitern.
equilibrium3.py
from ecell4 import *
with reaction_rules():
A + B == C | (0.01, 0.3)
m = get_model()
w = meso.MesoscopicWorld(Real3(1, 1, 1), Integer3(1, 1, 1)) # XXX: Point2
w.bind_to(m) # XXX: Point1
w.add_molecules(Species('C'), 60)
sim = meso.MesoscopicSimulator(w) # XXX: Point1
obs = FixedIntervalNumberObserver(0.1, ('A', 'B', 'C'))
sim.run(10, obs)
viz.plot_number_observer(obs)
Es ist fast dasselbe, außer dass Mesoscopic World und Mesoscopic Simulator verwendet werden, aber schauen wir uns einige neue Elemente an.
Zunächst verknüpfen wir in 'ww.to (m)' ein Modell mit World. Wir mussten dies in der Anfängerausgabe nicht tun, aber vergessen Sie nicht, dass die Attribute von Arten in räumlichen Techniken sehr wichtig sind. Stattdessen reicht es aus, dem Mesoscopic Simulator nur World zu geben.
Der nächste wichtige Unterschied besteht darin, dass wir bei der Erstellung der mesoskopischen Welt "Integeger3 (1, 1, 1)" als zweites Argument angeben. [^ 2] Dies wurde in ODE World oder Gillespie World nicht gesehen. Bevor wir erklären, was dies ist, ändern wir es und sehen, was passiert.
[^ 2]: 'Integer3' ist der gleiche 3D-Vektor wie 'Real3', aber das Element ist keine reelle Zahl Real, sondern ein ganzzahliger Wert Integer.
w = meso.MesoscopicWorld(Real3(1, 1, 1), Integer3(4, 4, 4)) # <- Integer3(1, 1, 1)
Ich denke, Sie haben ein anderes Ergebnis erzielt als zuvor. Je höher die Zahl, desto deutlicher der Unterschied.
Tatsächlich repräsentiert dieses zweite Argument die Anzahl der Unterteilungen im Raum. Meso ist ungefähr gleich wie Gillespie, aber während Gillespie die Berechnungen in einem einheitlichen geschlossenen Raum durchführt, Meso Teilt den Raum in kleine Quadrate, die als Subvolumes bezeichnet werden, so dass jedes Subvolume unterschiedliche molekulare Konzentrationen aufweisen kann. Daher wird im obigen Beispiel ein Würfel mit einer Seite in 64 Würfel mit einer Seite von 0,25 unterteilt.
Da zu Beginn 60 molekulare Spezies C eingegeben werden, ist höchstens 1 C in einem Subvolumen enthalten. Zu diesem Zeitpunkt ist die Konzentration von C dieselbe wie im ersten Fall, jedoch A nach der Trennung. Die Konzentration eines gerade abgetrennten A-Moleküls beträgt ursprünglich $ 1/1 = 1 $, im Subvolume jedoch 1 $ / (1/64) = 64 $. Mit dieser Anzahl von Teilungen wird der bisherige Effekt nicht erzielt, aber infolgedessen wird die Bindungsreaktion schneller sein, so dass die Menge des Bindungszustands C wie oben gezeigt zunimmt. Lassen Sie uns dies aus einer anderen Perspektive erklären. In meso Unmittelbar nach der Trennung befinden sich A und B immer im selben Subvolumen. Andererseits nimmt Gillespie ein gut gerührtes und gleichmäßiges System an, sodass A- und B-Moleküle unmittelbar nach der Trennung in einem der 64 Subvolumina enthalten sind. Zu diesem Zeitpunkt beträgt die Wahrscheinlichkeit, dass A und B in demselben Subvolumen enthalten sind, $ 1/64 $. Daher ist auch die Wahrscheinlichkeit einer Kombination hoch.
Gibt den Diffusionskoeffizienten des Moleküls an
Der Grund für diesen Unterschied ist, dass wir durch die Verwendung von Meso Platz bekommen haben, aber die Diffusionsbewegung des Moleküls noch nicht berücksichtigt haben.
Verwenden Sie dazu den Attributwert 'D' der Molekülspezies Species wie in der Anfängerausgabe. Als wir das letzte Mal add_species_attribute verwendet haben, aber wie in ReactionRule, ist hier die spezielle Notation für E-Cell 4. Ist verfügbar.
with species_attributes():
A | {'D': '1'}
B | {'D': '1'}
C | {'D': '1'}
# A | B | C | {'D': '1'} # means the same as above
Durch Verwendung von'species_attributes 'können Attribute mit dem Diktattyp [^ 3] nach dem Trennzeichen' | 'gesetzt werden. Hier wird' D ', was bedeutet, dass der Diffusionskoeffizient auf 1 gesetzt wird.
Fügen wir dies nun dem obigen Modell hinzu und berechnen es neu. Ich denke, es wird länger dauern als beim letzten Mal, aber ich denke, wir haben wieder ähnliche Ergebnisse wie Gilespie.
Warum löst die Diffusionsbewegung des Moleküls das vorherige Problem?
Betrachten Sie die freie Diffusionsbewegung eines Moleküls mit einem Diffusionskoeffizienten von $ D $ im dreidimensionalen Raum. Die Einheit des Diffusionskoeffizienten ist $ \ mathrm {\ mu m ^ 2 / s $ oder $ \ mathrm {nm} ^ 2 / \ Es ist das Quadrat der Länge geteilt durch die Zeit wie mu s $. Zu diesem Zeitpunkt ist bekannt, dass das Quadrat der Entfernung von der ersten Position zur Position nach der Zeit $ t $ Sekunden durchschnittlich $ 6Dt $ entspricht. Ist gewesen.
Umgekehrt beträgt die durchschnittliche Zeitskala für ein Molekül, um in einem Raum mit einer Längenskala von $ l $ zu bleiben, ungefähr $ l ^ 2 / 6D $. Im vorherigen Beispiel beträgt eine Seite des Subvolumens 0,25. Da der Diffusionskoeffizient 1 beträgt, beträgt er ungefähr 0,01 Sekunden. Selbst wenn sich A- und B-Moleküle im selben Subvolumen befinden, dauert die Reaktion ungefähr 1,5 Sekunden. Selbst wenn sie im selben Subvolumen getrennt sind, ist sie nahezu gleich. Im Fall von diffundiert es vor der Reaktion und bewegt sich zu einem anderen Subvolumen. Je kleiner das $ l $ ist, desto kleiner ist das Subvolumenvolumen $ l ^ 3 $ und desto höher ist daher die Reaktionsgeschwindigkeit nach der Trennung. Andererseits dauert es aufgrund der Diffusion weniger Zeit, um aus dem Raum herauszukommen.
Molekulare Lokalisation
Bisher wurde die Funktion'add_molecules 'verwendet, um der Welt Moleküle hinzuzufügen, wie bei Ode und Gillespie. Andererseits ermöglicht Mesoscopic World die Anordnung von Molekülen gemäß ihrer räumlichen Darstellung.
w = meso.MesoscopicWorld(Real3(1, 1, 1), Integer3(3, 3, 3))
w.add_molecules(Species('A'), 120)
w.add_molecules(Species('B'), 120, Integer3(1, 1, 1))
In Mesoscopic World'add_molecules 'können Sie angeben, in welches Subvolume das Molekül eingefügt werden soll, indem Sie nur eine Ganzzahl3 als drittes Argument angeben.
Im obigen Beispiel sind die Molekülspezies A im gesamten Raum verstreut, aber B befindet sich direkt im zentralen Subvolumen [^ 4]. Um dies zu bestätigen, verwenden Sie 'num_molecules'.
[^ 4]: Da die Koordinaten bei 0 beginnen, wird Integer3 (1, 1, 1) das zweite Subvolumen in jeder Achsenrichtung.
print(w.num_molecules(Species('B'))) # will print 120
print(w.num_molecules(Species('B'), Integer3(0, 0, 0))) # will print 0
print(w.num_molecules(Species('B'), Integer3(1, 1, 1))) # will print 120
Wenn Sie das IPython-Notizbuch verwenden, können Sie die Lokalisierung des Moleküls mithilfe der ecell4.viz-Bibliothek [^ 5] intuitiver erfassen.
viz.plot_world(w, radius=0.01)
'viz.plot_world' zeigt interaktiv die Anordnung der Partikel auf dem IPython-Notizbuch an, indem World [^ 6] angegeben wird. Legen Sie den Radius als Argument fest, um die Größe des Moleküls anzugeben.
[^ 6]: Das Molekül wird als Partikel angezeigt, aber in der mesoskopischen Welt gibt es keine bestimmte Position, an der es tatsächlich im Subvolumen vorhanden ist, sodass es jedes Mal zufällig im Subvolumen platziert wird.
Da wir die anfängliche Anordnung des Moleküls lokalisieren konnten, führen wir tatsächlich eine darauf basierende Simulation durch. Im obigen Beispiel wurde der Diffusionskoeffizient auf 1 und eine Seite der gesamten Welt auf 1 gesetzt, was ausreichend ist. 10 Sekunden reichen aus, um sich zu rühren. Rufen Sie nach der Berechnung erneut 'viz.plot_world' auf, um sicherzustellen, dass es tatsächlich bewegt wird.
Anfängliche Anordnung und Reaktion von Molekülen
Lassen Sie uns anhand eines extremen Beispiels sehen, wie sich die Lokalisierung des Moleküls auf die Reaktion auswirkt.
from ecell4 import *
with species_attributes():
A | B | C | {'D': '1'}
with reaction_rules():
A + B > C | 0.01
m = get_model()
w = meso.MesoscopicWorld(Real3(10, 1, 1), Integer3(10, 1, 1))
w.bind_to(m)
Das Modell enthält nur einfache Bindungen. World ist ein rechteckiger Körper, der nur in Richtung der X-Achse lang ist. Danach werden die Moleküle voreingenommen angeordnet.
w.add_molecules(Species('A'), 1200, Integer3(2, 0, 0))
w.add_molecules(Species('B'), 1200, Integer3(7, 0, 0))
viz.plot_world(w, radius=0.025)
Abgesehen davon gibt es einen Grund, hier nicht 'Ganzzahl3 (0, 0, 0)' und 'Ganzzahl3 (9, 0, 0)' zu verwenden. Grundsätzlich sind alle Welten in E-Zelle 4 periodische Randbedingungen [^ Periodische Grenze]. Daher sind die beiden vorherigen Subvolumes tatsächlich nebeneinander.
[^ PeriodicBoundary]: Periodische Randbedingung. Einfach ausgedrückt, wenn das Molekül auf eine bestimmte Achse fokussiert und das rechte Ende überschreitet, kommt das linke Ende heraus und umgekehrt.
Wenn die erwartete Anordnung erreicht ist, berechnen wir mit dem Mesoskopischen Simulator.
sim = meso.MesoscopicSimulator(w)
obs1 = NumberObserver(('A', 'B', 'C')) # XXX: saves the numbers after every steps
sim.run(5, obs1)
viz.plot_number_observer(obs1)
Ist die Reaktion ordnungsgemäß verlaufen? Um das Ergebnis zu bestätigen, versuchen Sie, die Anordnung gleichmäßig mit Meso zu verteilen, oder lassen Sie das Gillespie-Modul dasselbe Modell berechnen.
Wenn die durchgezogene Linie vorgespannt ist und keine gepunktete Linie vorhanden ist. Wenn eine Verzerrung vorliegt, tritt eine deutliche Verzögerung der Reaktion auf. Wenn Sie genauer hinschauen, stellen Sie möglicherweise fest, dass sich die Form der Zeitreihe auch zwischen der durchgezogenen Linie und der gepunkteten Linie unterscheidet. Nein. Dies liegt wahrscheinlich daran, dass es einige Zeit dauert, bis sich A- und B-Moleküle durch Diffusion treffen. Tatsächlich beträgt die Zeitskala, die erforderlich ist, um die Entfernung zwischen A und B (ungefähr 4) aufgrund der anfänglichen Anordnung zurückzulegen, $ 4 ^ 2/2. (D_ \ mathrm {A} + D_ \ mathrm {B}) = 4 $ Sekunden.
Einzelmolekül-Partikelgrößensimulation durch Spatiocyte
Wir haben ein Beispiel für die räumliche Darstellung in E-Cell 4 unter Verwendung der räumlichen Gillespie-Methode erläutert. Als nächstes betrachten wir die Berechnung der Molekülkorngröße 1, eine räumliche Darstellung mit höherer Auflösung.
Mikroskopische Gittermethode Spatiocyte
Bei der räumlichen Gillespie-Methode wird der Raum, der normalerweise von der Gillespie-Methode behandelt wird, in kleinere Räume unterteilt, und der Raum wird unter Berücksichtigung der Diffusion von Molekülen zwischen den Subvolumina einbezogen. Die Anzahl der Moleküle im Subvolumen ist jedoch immer noch groß. Es wird behandelt und die Position jedes Moleküls ist ungewiss. Mit anderen Worten, die räumliche Auflösung dieser Berechnungsmethode beträgt auf jeder Seite des Subvolumens $ l $. Um diese Auflösung zu erhöhen, fügen Sie $ l $ immer mehr hinzu. Es kann kleiner gemacht werden, aber bei dieser Methode muss $ l $ groß genug (mindestens 10-mal) in Bezug auf den Durchmesser des Moleküls $ R $ sein.
Wie können wir dann die räumliche Auflösung auf die Größe des Moleküls erhöhen? Die Antwort ist nicht die Anzahl der Moleküle, sondern die Simulation der Partikelgröße eines Moleküls, die die Reaktionsdiffusionsbewegung jedes Moleküls mit Position direkt ausdrückt. E-Cell 4 kann mehrere Einzelmolekül-Korngrößenberechnungstechniken verwenden, aber hier werden wir anhand der mikroskopischen Gittermethode Spatiocyte als Beispiel [^ Arjunan10] [^ Spatiocyte] erklären.
Spatiocyte behandelt jedes Molekül als eine reaktive starre Kugel und bewegt es in dem Raum, der durch das hexagonal dichteste Gitter [^ HCP] dargestellt wird. Im Gegensatz zur räumlichen Gillespie-Methode erhält jedes Molekül eine Identifikationsnummer. Anhand der Partikelgröße des Moleküls lässt sich erkennen, wo sich das Molekül befindet. In Bezug auf die Zeitskala ist die Zeitskala durch Diffusion proportional zum Quadrat der Entfernung, sodass sie 100-mal feiner ist als die räumliche Gillespie-Methode. Es wird die Zeitbreite sein.
[^ HCP]: [http://ja.wikipedia.org/wiki/Rokukata nächstgelegene Füllstruktur](http://ja.wikipedia.org/wiki/%E5%85%AD%E6%96%B9% E6% 9C% 80% E5% AF% 86% E5% 85% 85% E5% A1% AB% E6% A7% 8B% E9% 80% A0)
Wie auch immer, verwenden wir Spatiocyte. Die Grundlagen sind jedoch genau die gleichen.
__Hinweis: In den folgenden Beispielen wird das Namensgitter für Spatiocyte verwendet, es wird jedoch in einer zukünftigen Änderung in den Modulnamen Spatiocyte umbenannt.
equilibrium4.py
from ecell4 import *
with species_attributes():
A | B | C | {'D': '1'}
with reaction_rules():
A + B == C | (0.01, 0.3)
m = get_model()
w = lattice.LatticeWorld(Real3(1, 1, 1), 0.005) # The second argument is 'voxel_radius'.
w.bind_to(m)
w.add_molecules(Species('C'), 60)
sim = lattice.LatticeSimulator(w)
obs = FixedIntervalNumberObserver(0.1, ('A', 'B', 'C'))
sim.run(10, obs)
Der Unterschied ist das zweite Argument für LatticeWorld, das als Voxel-Radius bezeichnet wird. Spatiocyte unterteilt den Raum in Molekülgrößen, um die Position des Moleküls zu bestimmen, aber die kleinste Einheit dieses Raums ist Voxel. In den meisten Fällen kann dies die Größe des Moleküls sein. In diesem Beispiel wird der Radius des Moleküls im Raum von 1 $ \ mathrm {\ mu m} $ auf jeder Seite angegeben. 5 $ \ mathrm {nm} $.
Die Berechnungszeit wird länger als zuvor sein, aber die Ergebnisse stimmen natürlich mit Ode und Gillespie überein.
Diffusionsbewegung eines Moleküls
Lassen Sie uns tatsächlich ein Molekül behandeln, um die Berechnung der Partikelgröße eines Moleküls zu bestätigen.
from ecell4 import *
with species_attributes():
A | {'D': '1'}
m = get_model()
w = lattice.LatticeWorld(Real3(1, 1, 1), 0.005)
w.bind_to(m)
(pid, p), suc = w.new_particle(Species('A'), Real3(0.5, 0.5, 0.5))
Die Funktion 'new_particle' platziert ein Molekül in LatticeWorld an der angegebenen Stelle. Gleichzeitig gibt es die Identifikationsnummer' pid 'des erstellten Moleküls, das Tupel der Partikelinformation' p 'und'uc' zurück, ob es erfolgreich platziert wurde oder nicht. Wenn die Moleküle bereits an dieser Stelle platziert wurden, können sie nicht übereinander platziert werden, so dass suc 'falsch wird und versagt. Partikel' p 'enthält die Position, die Molekülspezies, den Radius und den Diffusionskoeffizienten des Moleküls. Von nun an kann das Partikel "p" extrahiert werden und die Information kann unter Verwendung der Identifikationsnummer "pid" erhalten werden.
Lass es uns überprüfen.
pid, p = w.get_particle(pid)
print(p.species().serial()) # will print: A
print(p.radius(), p.D()) # will print: (0.005, 1.0)
print(tuple(p.position())) # will print: (0.49806291436591293, 0.49652123150307814, 0.5)
Die Funktion'get_particle' empfängt die Identifikationsnummer und gibt die Identifikationsnummer und das Partikel zurück. Natürlich ist die Identifikationsnummer dieselbe wie die angegebene. Aus dem Partikel kann die Position durch die Funktion 'position' als Real3 abgerufen werden. Es ist schwierig zu lesen, wie sie ist. Es wird angezeigt, nachdem es in Tupel konvertiert wurde.
Die zurückgegebene Position liegt geringfügig von der angegebenen Position entfernt, da Spatiozyten nur Moleküle auf einem Gitter platzieren können. Lattice World platziert die Moleküle im nächstgelegenen Gitter von der angegebenen Position Wenn Sie die .plot_world'-Funktion verwenden, können Sie sehen, dass nur ein Molekül sicher zentriert ist.
Observer kann auch verwendet werden, um die Diffusion dieses Moleküls zu verfolgen.
sim = lattice.LatticeSimulator(w)
obs = FixedIntervalTrajectoryObserver(0.002, (pid,))
sim.run(1, obs)
viz.plot_trajectory(obs)
Hier wird es mit der Funktion 'viz.plot_trajectory' angezeigt, aber es ist auch möglich, die Trajektorie als Liste von Real3 durch die Funktion 'data' zu extrahieren, wie im Fall von NumberObserver.
print(len(obs.data())) # should return 1
print(len(obs.data()[0])) # should return 501
Die Datenfunktion gibt eine verschachtelte Liste zurück. Der Index in der ersten Liste ist der Index des Partikels. Da hier nur ein Partikel erzeugt wird, beträgt seine Indexlänge 1. Der Index in der nächsten Liste ist die Real3-Liste, die den Partikeln in der ersten Liste entspricht. Hier wurde die Simulation 1 Sekunde lang ausgeführt und die Flugbahn in Schritten von 0,002 Sekunden aufgezeichnet. Daher werden 501 Real3s zurückgegeben, einschließlich des Anfangswertes von einem Real3 und 1 / 0,002 = 500 Real3s.
Obwohl die Identifikationsnummer des durch die Funktion "add_molecules" hinzugefügten Moleküls unbekannt ist, können Partikel durch "list_particles" kollektiv aus den Molekülspezies extrahiert werden.
w.add_molecules(Species('A'), 5)
particles = w.list_particles(Species('A'))
for pid, p in particles:
print(p.species().serial(), tuple(p.position()))
Die Funktion 'list_particles' kann sowohl in anderen Welten als auch in der Funktion'add_molecules 'verwendet werden. [^ List_particles] Es ist wichtig, sich daran zu erinnern. Abgesehen davon ist die formale Funktion, die ein Molekül in Spatiocyte handhabt, ursprünglich' Es ist list_voxels ', und die Koordinaten werden durch einen ganzzahligen Wert anstelle von Real3 dargestellt.
[^ list_particles]: Es ist jedoch in ODE World nicht verfügbar.
Diffusionskoeffizient und Sekundärreaktion
Die Bindungsreaktion des Modells, mit dem wir uns bisher befasst haben, wird auch als Sekundärreaktion bezeichnet. Betrachten wir die Beziehung zwischen dieser Sekundärreaktion und dem Diffusionskoeffizienten in Spatiocyte.
from ecell4 import *
with species_attributes():
A | B | C | {'D': '1'}
with reaction_rules():
A + B > C | 1.0
m = get_model()
w = lattice.LatticeWorld(Real3(2, 1, 1), 0.005)
w.bind_to(m)
w.add_molecules(Species('A'), 120)
w.add_molecules(Species('B'), 120)
obs = FixedIntervalNumberObserver(0.005, ('A', 'B', 'C'))
sim = lattice.LatticeSimulator(w)
sim.run(1.0, obs)
%matplotlib inline
odew = ode.ODEWorld(Real3(2, 1, 1))
odew.bind_to(m)
odew.add_molecules(Species('A'), 120)
odew.add_molecules(Species('B'), 120)
odeobs = FixedIntervalNumberObserver(0.005, ('A', 'B', 'C'))
odesim = ode.ODESimulator(odew)
odesim.run(1.0, odeobs)
viz.plot_number_observer(obs, "-", odeobs, "--")
Dieses Mal verwenden wir eine höhere Geschwindigkeitskonstante als zuvor, aber wir können das gleiche Ergebnis wie zuvor erhalten. Wenn wir es jedoch im Detail mit dem Ergebnis der normalen Differentialgleichung vergleichen, können wir sehen, dass das Ergebnis unterschiedlich ist (die durchgezogene Linie in der folgenden Abbildung ist Spatiocyte). , Die gepunktete Linie ist auf Ode zurückzuführen).
Stimmt das nicht? Wenn die Geschwindigkeitskonstante auf einen viel höheren Wert eingestellt wird, stoppt die Reaktionsgeschwindigkeit endlos mit der Ode, aber mit Spatiocyte ist sie kaum schneller.
Dies ergibt sich aus dem Unterschied in der Definition der Reaktionsgeschwindigkeitskonstante zwischen dem Löser und der Einzelmolekül-Partikelgrößenberechnungstechnik. Die erstere wird als makroskopische oder effektive effektive Reaktionsgeschwindigkeitskonstante bezeichnet, während die letztere. Wird als mikroskopische oder intrinsische intrinsische Reaktionsgeschwindigkeitskonstante bezeichnet.
Die makroskopische Reaktionsgeschwindigkeitskonstante repräsentiert die Geschwindigkeit, mit der die Reaktion stattfindet, wenn die Moleküle gemischt werden, während die mikroskopische Reaktionsgeschwindigkeitskonstante den Grad der Reaktivität darstellt, wenn die Moleküle tatsächlich kollidieren. Aus mikroskopischer Sicht muss es daher zuerst kollidieren, bevor es reagiert. In Spatiocyte ist diese mikroskopische Geschwindigkeitskonstante, egal wie schnell sie ist, größer als die durch Diffusion verursachte Kollisionsgeschwindigkeit. Diese Situation wird als Diffusionsratenbestimmung bezeichnet. Dies ist sehr ähnlich zu der Zeit, die die vorgespannten Moleküle brauchten, um bei der räumlichen Gillespie-Methode miteinander zu reagieren.
Die folgende Beziehung ist im dreidimensionalen Raum zwischen dieser makroskopischen Geschwindigkeitskonstante $ k_ \ mathrm {on} $ und der mikroskopischen Geschwindigkeitskonstante $ k_a $ [^ Collins Kimball] bekannt.
\frac{1}{k_\mathrm{on}}=\frac{1}{k_a}+\frac{1}{4\pi RD_\mathrm{tot}}
Dabei ist $ R $ die Summe der Radien der beiden kollidierenden Moleküle und $ D_ \ mathrm {tot} $ die Summe der Diffusionskoeffizienten. Im obigen Beispiel ist der zweite Term auf der rechten Seite $ k_D = 4 \ pi RD_ \ Da mathrm {tot} $ ungefähr 0,25 beträgt, beträgt die makroskopische Reaktionsgeschwindigkeitskonstante ungefähr 0,2, wenn sie mit der mikroskopischen Reaktionsgeschwindigkeitskonstante 1,0 & agr; kombiniert wird.
[^ alpha]: In Spatiocyte muss die Konstante der sekundären Reaktionsgeschwindigkeit jedoch langsamer als $ 3 \ sqrt {2} RD $ sein, und die Kollisionsratenkonstante $ k_D $ beträgt ebenfalls $ 3 \ sqrt, sofern keine besonderen Einstellungen vorgenommen werden. {2} RD $.
Im Gegensatz zu den normalen Differentialgleichungen und der Gillespie-Methode, die ein derart gut gerührtes System voraussetzen (dh der Diffusionskoeffizient ist unendlich), kann die Berechnung der Partikelgröße eines Moleküls die Diffusion und Reaktion des Moleküls richtig trennen. Wie aus der Gleichung ersichtlich ist, stimmt die mikroskopische Geschwindigkeitskonstante, wenn sie in Bezug auf $ k_D $ ausreichend klein ist, nahezu mit der makroskopischen Geschwindigkeitskonstante überein (Reaktionsgeschwindigkeitsbestimmung).
Behandlung von Strukturen in Spatiocyte
Abschließend werde ich erklären, wie mit Strukturen wie Zellmembranen mit Spatiocyte umgegangen wird. In E-Cell 4 befindet sich die Struktur noch in der Entwicklungsphase, aber Spatiocyte kann ihre Funktion in gewissem Umfang nutzen.
Derzeit sind die Formen, die verarbeitet werden können, begrenzt. Betrachten wir jedoch zunächst eine Kugel als Beispiel. Um die Bewegung von Molekülen auf einen kugelförmigen Raum zu beschränken, erstellen Sie zunächst eine kugelförmige Struktur.
w = lattice.LatticeWorld(Real3(1, 1, 1), 0.005)
sph = Sphere(Real3(0.5, 0.5, 0.5), 0.45)
print(w.add_structure(Species('C'), sph)) # will print 539805
viz.plot_world(w)
'Kugel' bedeutet eine Kugel mit dem ersten Argument als Zentrum und dem zweiten Argument als Radius. Ich gab ihr die Spezies 'C' und fügte sie der Welt hinzu.
In Spatiocyte wird die Struktur dargestellt, indem der Bereich im Raum mit einem speziellen Voxel gefüllt wird. Im obigen Beispiel sind alle Voxel in der angegebenen Sphäre von der Molekülspezies C.'viz besetzt. Wenn Sie es mit plot_world 'anzeigen, können Sie sehen, dass es tatsächlich sphärisch verteilt ist. Da die Anzahl der Moleküle jedoch zu groß ist, um etwa 540.000 anzuzeigen, wird nur ein Teil geplottet, aber tatsächlich Ist komplett gefüllt.
Erstellen Sie nach dem Erstellen der Kugel ein Molekül, das sich nur in dieser Kugel bewegt. Geben Sie dazu das Attribut "Position" an.
with species_attributes():
A | {'D': '1', 'location': 'C'}
m = get_model()
Dann'add_molecules 'die Moleküle wie gewohnt.
w.bind_to(m)
w.add_molecules(Species('A'), 120)
viz.plot_world(w, species_list=('A',)) # visualize A-molecules only
Nachdem wir nun angegeben haben, dass der molekulare Typ A auf der Struktur des molekularen Typs C existiert, wird sogar 'add_molecules' tatsächlich nur in der Kugel platziert. Es sollte beachtet werden, dass natürlich die Struktur vor dem Hinzufügen von A vorhanden ist. Müssen gemacht werden.
Lassen Sie uns überprüfen, ob die Bewegung dieser molekularen Spezies A in der Sphäre wirklich eingeschränkt ist. Auch hier kann 'FixedIntervalTrajectoryObserver' verwendet werden.
pid_list = [pid for pid, p in w.list_particles(Species('A'))[: 10]]
obs = FixedIntervalTrajectoryObserver(1e-3, pid_list)
sim = lattice.LatticeSimulator(w)
sim.run(1, obs)
viz.plot_trajectory(obs)
'pid_list' ist eine Liste von 10 geeigneten Identifikationsnummern von 60 A-Molekülen. Farbig die Diffusionsbahn dieser 10 Moleküle. Sicherlich ist die Bewegung in der Kugel eingeschränkt. Ich verstehe.
Struktur und Reaktion
Lassen Sie mich abschließend die Migration von Molekülen zwischen Strukturen erklären. Erstens gehören molekulare Spezies ohne das Attribut "Ort" zu keiner Struktur. Geben Sie im vorherigen Beispiel "Ort" an. Andernfalls wird es außerhalb der Kugel platziert.
Nehmen wir als Beispiel eine Ebene und probieren Sie es aus. Um eine Ebene zu erstellen, können Sie sie zunächst wie folgt erstellen, indem Sie drei Real3s verwenden, den Ursprung 'Ursprung' und die zwei Achsenvektoren 'Einheit0' und 'Einheit1'.
ps = PlanarSurface(origin, unit0, unit1)
Unter Verwendung eines ebenen molekularen Typs A und eines normalen molekularen Typs B,
from ecell4 import *
with species_attributes():
A | {'D': '0.1', 'location': 'M'}
B | {'D': '1'}
m = get_model()
w = lattice.LatticeWorld(Real3(1, 1, 1))
w.bind_to(m)
origin = Real3(0, 0, 0.5)
unit0 = Real3(1, 0, 0)
unit1 = Real3(0, 1, 0)
w.add_structure(
Species('M'), PlanarSurface(origin, unit0, unit1)) # Create a structure first
w.add_molecules(Species('B'), 480) # Throw-in B-molecules
viz.plot_world(w, species_list=('A', 'B'))
Es ist ein wenig schwer zu verstehen, aber B-Moleküle sind nur in der Ebene angeordnet. Wie sollen wir dann die Reaktion ausdrücken, bei der diese B-Moleküle in der Ebene M adsorbiert werden und zu A-Molekülen werden?
with reaction_rules():
B + M == A | (1e-3, 1.5)
Wie Sie sehen können, wird das B-Molekül, wenn es mit der Struktur M kollidiert, adsorbiert und wird zum A-Molekül. Im Gegenteil, das A-Molekül wird zur ursprünglichen Ebene und das B-Molekül, wenn es getrennt wird.
Zu diesem Zeitpunkt können Sie wie gewohnt mit Simulator rechnen.
sim = lattice.LatticeSimulator(w)
obs = NumberObserver(('A', 'B'))
sim.run(2, obs)
viz.plot_number_observer(obs)
viz.plot_world(w, species_list=('A', 'B'))
Übrigens kann die molekulare Spezies der Struktur weggelassen werden, indem sie von der Struktur getrennt wird, aber sie kann nicht durch Bindung weggelassen werden. Es ist unmöglich, A aus B zu machen, wenn es kein M gibt, um das das Produkt-A-Molekül gelegt werden sollte. Das ist der Grund. Im Gegenteil, in diesem Fall gibt es kein Problem, da das A-Molekül das B-Molekül irgendwo auf der Ebene M platzieren kann. Außerdem ist im Fall der Trennung die Struktur primär oder nicht weggelassen. Es ist immer noch eine Reaktion, aber bei der Bindung wird die Sekundärreaktion zu einer Primärreaktion, und auch die Bedeutung der Geschwindigkeitskonstante ändert sich.
with reaction_rules():
B + M > A | 1e-3
A > B | 3.0 # means the same as A > B + M
Übungen
-
Übungen zur räumlichen Gillespie-Methode: Die räumliche Gillespie-Methode wurde verwendet, um die Zeitskala der Diffusion für ein Subvolumen mit einer Seite von $ l $ zu berechnen. Andererseits wurde auch die Zeitskala erläutert, die für die Moleküle unmittelbar nach der Trennung erforderlich ist, um sie wieder anzubringen. Betrachten Sie die Abhängigkeit dieser beiden Zeitskalen von $ l $ und klären Sie die Bedingungen, die $ l $ erfüllen sollten. Es ist jedoch bekannt, dass die räumliche Gillespie-Methode die folgenden Bedingungen für $ l $ erfüllen sollte. Ja.
R^2 \ll l^2 \ll 6D\tau_\mathrm{min}
Dabei ist $ R $ der Durchmesser des Moleküls (genauer gesagt der Reaktionsradius) und $ \ tau_ \ mathrm {min} $ die kürzeste Lebensdauer für jede Reaktion im Molekül [^ Elf04].
[^Elf04]: Elf J. and Ehrenberg M., "Spontaneous separation of bi-stable biochemical systems into spatial domains of opposite phases", Syst. Biol., 1(2), 230-236 (2004)
1. Übung zur anfänglichen Anordnung von Molekülen: Ändern Sie das obige Beispiel, um zu untersuchen, wie sich die anfängliche Anordnung von AB-Molekülen, Reaktionsgeschwindigkeitskonstanten, Diffusionskoeffizienten und anderen Parametern auf die Reaktion auswirkt. Ist es möglich, ein Modell zu erstellen, das die Zeitreihen wie im obigen Beispiel mit der einheitlichen Gillespie-Methode reproduziert?
1. Übung zur molekularen Diffusion: Für die freie Diffusion in drei Dimensionen ist die mittlere quadratische Verschiebung der durchschnittliche quadratische Abstand zur Position des Moleküls mit einem Diffusionskoeffizienten von $ D $ nach $ t $ Sekunden.
```math
\left<L^2\right>=6Dt
Ist bekannt, dass dies der Fall ist. Ist dies tatsächlich bei Spatiocyte der Fall?
-
Übungen zu makroskopischen und mikroskopischen Geschwindigkeitskonstanten: Sowohl makroskopische als auch mikroskopische Geschwindigkeitskonstanten haben die gleiche Einheit von $ \ mathrm {nM} ^ {-1} s ^ {-1} $ und $ \ mathrm {\ mu m} ^ 3s ^ {-1} $. Zwei Moleküle mit einer Diffusionsrate von 1 $ \ mathrm {\ mu m} ^ 2 / s $ und einem Durchmesser von 5 $ \ mathrm {nm} $ Was ist die Kollisionsratenkonstante $ k_D $ $ \ mathrm {nM} ^ {-1} s ^ {-1} $?
-
Übung zur Struktur: Überprüfen Sie mit einer flachen Struktur die Beziehung zwischen dem Diffusionskoeffizienten $ D $ in der Ebene und seinem durchschnittlichen quadratischen Abstand.
Recommended Posts