[PYTHON] E-Cell 4 Intermédiaire facile à utiliser

Préparation

Veuillez vous référer à l'article séparé pour l'installation.

• Essayez d'utiliser E-Cell 4 sous Windows 7 ou Mac OS X

Cet article est une édition intermédiaire. On suppose que vous comprenez l'édition pour débutants.

• E-Cell 4 facile à utiliser pour les débutants

De plus, cet article suppose l'utilisation d'IPython Notebook.

Comment utiliser

Espace dans E-Cell 4

Dans l'édition pour débutants, j'ai expliqué les trois éléments Model, World et Simulator dans E-Cell 4. J'ai également essayé un calcul simple en utilisant ode, qui est un solveur d'équation différentielle normal, et gilespie, qui est une méthode probabiliste.

J'ai créé un monde en donnant une taille lors de l'exécution de calculs avec ode et gilespie, mais quelle est la gestion de l'espace dans E-Cell 4?

from ecell4 import *

w1 = ode.ODEWorld(Real3(1, 1, 1))
w2 = gillespie.GillespieWorld(Real3(1, 1, 1))

Dans ode et gillespie, nous avons calculé dans un cube avec un côté comme dans l'exemple ci-dessus, mais en réalité seul ce volume compte.

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'

Cependant, comme 1 est après tout donné comme volume, le même résultat de calcul sera affiché.

Cela semble en fait raisonnable dans un système bien agité et spatialement uniforme, comme in vitro, mais les cellules ne sont clairement pas spatialement uniformes, pour tenir compte de la localisation de ces molécules. Nécessite des calculs biochimiques spatiaux.

Différentes expressions spatiales et techniques de calcul correspondantes peuvent être utilisées dans E-Cell 4. Ci-dessous, nous allons d'abord examiner la gestion de l'espace dans E-Cell 4 en utilisant la méthode spatiale de Gillespie, qui est l'une d'entre elles, à titre d'exemple.

Méthode Spatial Gillespie

Dans E-Cell 4, la méthode spatiale de Gillespie [^ 1] est incluse dans le module "méso", calculons d'abord en utilisant "run_simulation" comme nous l'avons vu dans l'édition débutant sans penser à rien.

[^ 1]: Plus correctement, il est plus proche de la méthode de sous-volume suivant. Pour plus de détails, reportez-vous au document.

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

Cela peut avoir pris un certain temps par rapport à ode et gillespie, mais vous obtiendrez presque le même résultat. C'est vrai, le module meso fait presque le même calcul que gilespie sans aucun paramètre spatial. Développons "run_simulation".

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)

C'est presque la même chose sauf qu'il utilise Mesoscopic World et Mesoscopic Simulator, mais jetons un coup d'œil à quelques nouveaux éléments.

Tout d'abord, dans "w.bind_to (m)", nous associons un modèle à World, ce que nous n'avons pas eu à faire dans l'édition débutant, mais n'oubliez pas que les attributs de Species sont très importants dans les techniques spatiales. Au lieu de cela, il suffit de donner au simulateur mésoscopique uniquement du monde.

La prochaine différence importante est que lors de la création du monde mésoscopique, nous donnons «Integer3 (1, 1, 1)» comme deuxième argument. [^ 2]. Cela n'a pas été vu dans ODE World ou Gillespie World. Avant d'expliquer ce que c'est, changeons-le et voyons ce qui se passe.

[^ 2]: "Integer3" est le même vecteur 3D que "Real3", mais l'élément n'est pas un nombre réel Real mais une valeur entière Integer.

w = meso.MesoscopicWorld(Real3(1, 1, 1), Integer3(4, 4, 4))  # <- Integer3(1, 1, 1)

Je pense que nous avons obtenu un résultat différent de celui d'avant: plus le nombre est élevé, plus la différence sera notable.

En fait, ce second argument représente le nombre de divisions dans l'espace. Meso est presque le même que gillespie, mais tandis que gillespie effectue les calculs dans un espace fermé uniforme, meso Divise l'espace en petits carrés appelés sous-volumes, permettant à chaque sous-volume d'avoir des concentrations moléculaires différentes. Par conséquent, dans l'exemple ci-dessus, un cube avec un côté est divisé en 64 cubes avec un côté de 0,25.

Étant donné que 60 espèces moléculaires C sont entrées au début, au plus 1 C est contenu dans un sous-volume.A ce moment, la concentration de C est la même que dans le premier cas, mais A après séparation. La concentration d'une molécule A qui vient d'être séparée est à l'origine $ 1/1 = 1 $, mais en Subvolume elle devient $ 1 / (1/64) = 64 $. Avec ce nombre de divisions, l'effet jusqu'à présent ne sera pas obtenu, mais en conséquence, la réaction de liaison sera plus rapide, de sorte que la quantité d'état de liaison C augmentera comme indiqué ci-dessus. Expliquons sous un autre angle. In meso Immédiatement après la séparation, A et B sont toujours dans le même sous-volume.En revanche, gillespie suppose un système bien agité et uniforme, de sorte que les molécules A et B immédiatement après la séparation sont incluses dans l'un des 64 sous-volumes. À ce stade, la probabilité que A et B soient inclus dans le même sous-volume est de 1/64 $. Par conséquent, la probabilité de combinaison est également élevée.

Donne le coefficient de diffusion de la molécule

La raison de cette différence est que nous avons obtenu de l'espace en utilisant le méso, mais nous n'avons pas encore considéré le mouvement de diffusion de la molécule.

Pour définir cela, utilisez l'attribut value'D 'de l'espèce moléculaire Species comme vous l'avez fait dans l'édition pour débutants. La dernière fois que nous avons utilisé add_species_attribute, mais comme nous l'avons fait dans ReactionRule, voici la notation spéciale pour E-Cell 4. Est disponible.

with species_attributes():
    A | {'D': '1'}
    B | {'D': '1'}
    C | {'D': '1'}

    # A | B | C | {'D': '1'}  # means the same as above

En utilisant 'species_attributes', vous pouvez définir des attributs avec le type de dict [^ 3] après le séparateur '|'. Ici, 'D' qui signifie que le coefficient de diffusion est mis à 1.

Maintenant, ajoutons ceci au modèle ci-dessus et recalculons. Je pense que le calcul prendra plus de temps que la dernière fois, mais je pense que nous avons à nouveau des résultats similaires à ceux de gilespie.

Alors pourquoi le mouvement de diffusion de la molécule résout-il le problème précédent?

Considérons le mouvement de diffusion libre d'une molécule avec un coefficient de diffusion de $ D $ dans un espace tridimensionnel. L'unité du coefficient de diffusion est $ \ mathrm {\ mu m} ^ 2 / s $ ou $ \ mathrm {nm} ^ 2 / . C'est le carré de la longueur divisé par le temps comme mu s $. A ce moment, on sait que le carré de la distance entre la première position et la position après le temps $ t $ secondes est égal à $ 6Dt $ en moyenne. A été.

Inversement, l'échelle de temps moyenne à laquelle une molécule reste dans un espace avec une échelle de longueur de $ l $ est d'environ $ l ^ 2 / 6D $. Dans l'exemple précédent, un côté du sous-volume est de 0,25. Puisque le coefficient de diffusion est de 1, il est d'environ 0,01 seconde. Par contre, même si les molécules A et B sont dans le même sous-volume, il faut environ 1,5 seconde pour réagir, donc même si elles sont séparées dans le même sous-volume, c'est presque la même chose. Dans le cas de, il diffuse avant la réaction et se déplace vers un autre sous-volume. En revanche, il faut moins de temps pour sortir de l'espace en raison de la diffusion.

Localisation moléculaire

Jusqu'à présent, la fonction "add_molecules" a été utilisée pour ajouter des molécules au Monde, comme avec ode et gillespie. D'autre part, Mesoscopic World permet l'arrangement des molécules selon leur représentation spatiale.

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))

Dans Mesoscopic World'add_molecules ', vous pouvez spécifier dans quel sous-volume placer la molécule en ne donnant qu'un seul Integer3 comme troisième argument.

Dans l'exemple ci-dessus, les espèces moléculaires A sont dispersées dans tout l'espace, mais B est situé juste dans le sous-volume central [^ 4]. Pour le confirmer, utilisez "num_molécules".

[^ 4]: Puisque les coordonnées commencent à 0, Integer3 (1, 1, 1) devient le deuxième sous-volume dans chaque direction d'axe.

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

De plus, si vous utilisez le notebook IPython, vous pouvez appréhender la localisation de la molécule de manière plus intuitive en utilisant la bibliothèque ecell4.viz [^ 5].

viz.plot_world(w, radius=0.01)

'viz.plot_world' affichera de manière interactive la disposition des particules sur le notebook IPython en donnant World [^ 6]. Définissez le rayon comme argument pour spécifier la taille de la molécule.

[^ 6]: La molécule est affichée sous forme de particule, mais dans Mesoscopic World, il n'y a pas de position spécifique où elle existe réellement dans le sous-volume, elle est donc placée au hasard dans le sous-volume à chaque fois.

Puisque nous avons pu localiser l'arrangement initial de la molécule, effectuons en fait une simulation basée sur cela. Dans l'exemple ci-dessus, le coefficient de diffusion a été fixé à 1 et un côté du monde entier a été mis à 1, donc c'est suffisant. 10 secondes suffiront pour remuer. Après le calcul, appelez à nouveau "viz.plot_world "pour vous assurer qu'il est bien agité.

Arrangement initial et réaction des molécules

Voyons comment la localisation de la molécule affecte la réaction avec un exemple extrême.

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)

Le modèle ne contient que des liaisons simples. Le monde est un corps rectangulaire qui n'est long que dans la direction de l'axe X. Ensuite, les molécules sont disposées de manière biaisée.

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)

En passant, il y a une raison de ne pas utiliser ici "Entier3 (0, 0, 0)" et "Entier3 (9, 0, 0)". Fondamentalement, tous les mondes de la cellule E 4 sont des conditions aux limites périodiques [^ Période de limite]. Par conséquent, les deux sous-volumes précédents sont en fait côte à côte.

[^ PeriodicBoundary]: condition aux limites périodiques. En termes simples, lorsque l'on se concentre sur un certain axe, lorsque la molécule dépasse l'extrémité droite, l'extrémité gauche sort, et vice versa.

Si l'arrangement attendu est obtenu, calculons à l'aide du simulateur mésoscopique.

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)

La réaction s'est-elle déroulée correctement? Pour confirmer le résultat, essayez de répartir l'arrangement uniformément avec meso, ou laissez le module gillespie calculer le même modèle.

Lorsque le trait plein est biaisé et qu'il n'y a pas de ligne pointillée. Lorsqu'il y a un biais, il y a un net retard dans la réaction. Si vous regardez de plus près, vous remarquerez peut-être que la forme de la série chronologique est également différente entre la ligne continue et la ligne pointillée. Non. C'est probablement parce que les molécules A et B mettent du temps à se rencontrer par diffusion. En fait, l'échelle de temps nécessaire pour parcourir la distance entre A et B (environ 4) en raison de la disposition initiale est de 4 $ ^ 2/2. (D_ \ mathrm {A} + D_ \ mathrm {B}) = 4 $ secondes.

Simulation de la taille des particules d'une molécule unique par Spatiocyte

Nous avons expliqué un exemple de représentation spatiale dans E-Cell 4 en utilisant la méthode spatiale de Gillespie, puis examinons le «calcul de la taille de grain moléculaire 1», qui est une représentation spatiale à plus haute résolution.

Méthode de réseau microscopique Spatiocyte

Dans la méthode spatiale de Gillespie, l'espace normalement manipulé par la méthode de Gillespie est divisé en espaces plus petits, et l'espace est introduit en considérant la diffusion de molécules entre les sous-volumes.Toutefois, le nombre de molécules dans le sous-volume est encore important. Elle est traitée et la position de chaque molécule est incertaine. Autrement dit, la résolution spatiale de cette méthode de calcul est égale à $ l $ de chaque côté du sous-volume. Pour augmenter cette résolution, $ l $ est de plus en plus ajouté. Vous pouvez le rendre plus petit, mais avec cette méthode, $ l $ doit être suffisamment grand (au moins 10 fois plus grand) que le diamètre de la molécule $ R $.

Alors, comment augmenter la résolution spatiale à la taille de la molécule? La réponse n'est pas le nombre de molécules, mais la simulation de taille de particule à une molécule qui exprime directement le mouvement de diffusion de réaction de chaque molécule avec position. Dans E-Cell 4, plusieurs techniques de calcul de la taille des particules à une seule molécule peuvent être utilisées, mais ici nous allons expliquer en utilisant la méthode du réseau microscopique Spatiocyte comme exemple [^ Arjunan10] [^ Spatiocyte].

Spatiocyte traite chaque molécule comme une sphère rigide réactive et la déplace dans l'espace représenté par le réseau hexagonal le plus proche [^ HCP] Contrairement à la méthode spatiale de Gillespie, chaque molécule reçoit un numéro d'identification. Il est possible de savoir où se trouve la molécule par la granulométrie de la molécule.En termes d'échelle de temps, l'échelle de temps due à la diffusion est proportionnelle au carré de la distance, elle est donc 100 fois plus fine que la méthode spatiale de Gillespie. Ce sera la largeur du temps.

[^ HCP]: [http://ja.wikipedia.org/wiki/Rokukata structure densément compactée](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)

Bref, utilisons Spatiocyte, mais les bases sont exactement les mêmes.

Remarque: Dans les exemples suivants, le nom treillis est utilisé pour Spatiocyte, mais il sera renommé en nom de module spatiocyte lors d'un changement ultérieur .

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)

La différence est le deuxième argument donné à LatticeWorld, qui est appelé le rayon de Voxel. Spatiocyte divise l'espace en tailles moléculaires pour déterminer la position de la molécule, mais la plus petite unité de cet espace est Voxel. Dans la plupart des cas, cela peut être la taille de la molécule.Dans cet exemple, le rayon de la molécule est donné dans l'espace de 1 $ \ mathrm {\ mu m} $ de chaque côté. 5 $ \ mathrm {nm} $.

Le temps de calcul sera plus long qu'auparavant, mais les résultats seront bien sûr cohérents avec ode et gillespie.

Mouvement de diffusion d'une molécule

Traitons en fait d'une molécule pour confirmer le calcul de la taille des particules d'une molécule.

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))

La fonction "nouvelle_particule" place une molécule dans LatticeWorld à l'emplacement spécifié. En même temps, elle renvoie le numéro d'identification "pid" de la molécule créée, le tuple de l'information de particule "p" et "uc "si elle a été placée avec succès. Si les molécules sont déjà placées à cet endroit, elles ne peuvent pas être placées les unes sur les autres, donc «u» devient Faux et échoue. La particule «p» contient la position, l'espèce moléculaire, le rayon et le coefficient de diffusion de la molécule. Désormais, la particule «p» peut être extraite et les informations peuvent être obtenues en utilisant le numéro d'identification «pid».

Regardons ça.

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)

La fonction "obtenir_particule" reçoit le numéro d'identification et renvoie le numéro d'identification et la particule. Naturellement, le numéro d'identification est le même que celui donné. A partir de la particule, la position peut être récupérée comme Réel3 par la fonction "position". Il est difficile à lire tel quel. Il s'affiche après avoir été converti en tuple.

La position renvoyée est légèrement en dehors de l'emplacement spécifié car les spatiocytes ne peuvent placer des molécules que sur une grille. Lattice World placera les molécules dans la grille la plus proche de l'emplacement donné. Si vous utilisez la fonction .plot_world ', vous pouvez voir qu'une seule molécule est centrée.

Observer peut également être utilisé pour suivre la diffusion de cette molécule.

sim = lattice.LatticeSimulator(w)
obs = FixedIntervalTrajectoryObserver(0.002, (pid,))
sim.run(1, obs)
viz.plot_trajectory(obs)

Ici, il est affiché à l'aide de la fonction'viz.plot_trajectory ', mais il est également possible d'extraire la trajectoire sous forme de liste de Real3 par la fonction'data' comme dans le cas de NumberObserver.

print(len(obs.data())) # should return 1
print(len(obs.data()[0])) # should return 501

La fonction de données renvoie une liste imbriquée. L'index de la première liste est l'indice de la particule. Puisqu'une seule particule est créée ici, sa longueur d'index est de 1. L'index de la liste suivante est la liste Real3 qui correspond aux particules de la première liste. Ici, la simulation a été exécutée pendant 1 seconde et la trajectoire a été enregistrée par incréments de 0,002 seconde. Par conséquent, 501 Real3 sont renvoyés, y compris la valeur initiale d'un Real3 et 1 / 0,002 = 500 Real3.

De plus, bien que le numéro d'identification de la molécule ajoutée par la fonction "add_molecules" soit inconnu, les particules peuvent être collectivement extraites de l'espèce moléculaire par "list_particles".

w.add_molecules(Species('A'), 5)

particles = w.list_particles(Species('A'))
for pid, p in particles:
    print(p.species().serial(), tuple(p.position()))

La fonction 'list_particles' peut être utilisée dans d'autres mondes ainsi que la fonction 'add_molecules' [^ list_particles] Il est important de se souvenir. En passant, la fonction formelle qui gère une molécule dans Spatiocyte est à l'origine ' Il s'agit de list_voxels 'et les coordonnées sont représentées par une valeur entière au lieu de Real3.

[^ list_particles]: Cependant, il n'est pas disponible dans ODE World.

Coefficient de diffusion et réaction secondaire

La réaction de liaison du modèle que nous avons traité jusqu'à présent est également appelée réaction secondaire.Examinons la relation entre cette réaction secondaire et le coefficient de diffusion chez 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, "--")

Cette fois, nous utilisons une constante de vitesse plus élevée qu'avant, mais nous pouvons obtenir le même résultat qu'auparavant. Cependant, si nous le comparons en détail avec le résultat de l'équation différentielle normale, nous pouvons voir que le résultat est différent (le trait plein dans la figure ci-dessous est Spatiocyte). , La ligne pointillée est due à ode).

En fait, lorsque la constante de vitesse est réglée sur une valeur beaucoup plus élevée, la vitesse de réaction s'arrête indéfiniment avec ode, mais elle n'est guère plus rapide avec Spatiocyte.

Cela est dû à la différence de définition de la constante de vitesse de réaction entre le solveur et la technique de calcul de la taille des particules à une seule molécule. La première est appelée constante de vitesse de réaction effective macroscopique ou effective, tandis que la seconde. Est appelée la constante de vitesse de réaction microscopique ou intrinsèque intrinsèque.

La constante de vitesse de réaction macroscopique représente la vitesse à laquelle la réaction se produit lorsque les molécules sont mélangées, tandis que la constante de vitesse de réaction microscopique représente le degré de réactivité lorsque les molécules entrent réellement en collision. Par conséquent, d'un point de vue microscopique, il doit d'abord entrer en collision avant de réagir. Cette situation s'appelle la détermination de la vitesse de diffusion, ce qui est très similaire au temps qu'il a fallu aux molécules biaisées pour réagir entre elles dans la méthode spatiale de Gillespie.

La relation suivante est connue dans l'espace tridimensionnel entre cette constante de vitesse macroscopique $ k_ \ mathrm {on} $ et la constante de vitesse microscopique $ k_a $ [^ Collins Kimball].

[^ CollinsKimball]: formule Smoluchowski – Collins – Kimball.

\frac{1}{k_\mathrm{on}}=\frac{1}{k_a}+\frac{1}{4\pi RD_\mathrm{tot}}

Où $ R $ est la somme des rayons des deux molécules en collision et $ D_ \ mathrm {tot} $ est la somme des coefficients de diffusion. Dans l'exemple ci-dessus, le deuxième terme du côté droit $ k_D = 4 \ pi RD_ \ Puisque mathrm {tot} $ est d'environ 0,25, la constante de vitesse de réaction macroscopique est d'environ 0,2 lorsqu'elle est combinée avec la constante de vitesse de réaction microscopique 1,0 [^ alpha].

[^ alpha]: Cependant, dans Spatiocyte, la constante de vitesse de réaction secondaire doit être plus lente que $ 3 \ sqrt {2} RD $, et la constante de taux de collision $ k_D $ est également $ 3 \ sqrt, sauf si des paramètres spéciaux sont définis. {2} RD $.

Contrairement aux équations différentielles normales et à la méthode de Gillespie, qui supposent un système aussi bien agité (c'est-à-dire que le coefficient de diffusion est infini), le calcul de la taille des particules à une molécule peut séparer correctement la diffusion et la réaction des molécules. Comme le montre l'équation de, si la constante de vitesse microscopique est suffisamment petite par rapport à $ k_D $, elle correspond presque à la constante de vitesse macroscopique (détermination de la vitesse de réaction).

Traitement des structures chez Spatiocyte

Enfin, j'expliquerai comment manipuler des structures telles que les membranes cellulaires en utilisant Spatiocyte.Dans E-Cell 4, la structure est encore au stade de développement, mais Spatiocyte peut utiliser sa fonction dans une certaine mesure.

Actuellement, les formes manipulables sont limitées, mais regardons d'abord une sphère à titre d'exemple.Pour limiter le mouvement des molécules à un espace sphérique, créez d'abord une structure sphérique.

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)

«Sphère» signifie une sphère avec le premier argument comme centre et le second argument comme rayon. Je lui ai donné une espèce appelée «C» et l'ai ajoutée à World.

Dans Spatiocyte, la structure est représentée en remplissant la zone de l'espace avec un Voxel spécial.Dans l'exemple ci-dessus, tous les Voxels de la sphère spécifiée sont occupés par l'espèce moléculaire C.'viz. Si vous l'affichez avec plot_world ', vous pouvez voir qu'il est en fait distribué sous une forme sphérique. Cependant, comme le nombre de molécules est trop grand pour être affiché à environ 540 000, seule une partie est tracée, mais en réalité Est complètement rempli.

Une fois la sphère créée, créez une molécule qui se déplace uniquement dans cette sphère. Pour ce faire, spécifiez l'attribut 'emplacement'.

with species_attributes():
    A | {'D': '1', 'location': 'C'}

m = get_model()

Puis "ajoutez_molécules" les molécules comme d'habitude.

w.bind_to(m)
w.add_molecules(Species('A'), 120)
viz.plot_world(w, species_list=('A',))  # visualize A-molecules only

Maintenant que nous avons précisé que le type moléculaire A existe sur la structure de type moléculaire C, même "add_molecules" est en fait placé uniquement dans la sphère. Il faut noter que, bien sûr, la structure avant d'ajouter A. Besoin d'être fait.

Vérifions si le mouvement de cette espèce moléculaire A est vraiment restreint dans la sphère. Là encore, 'FixedIntervalTrajectoryObserver' peut être utilisé.

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' est une liste de 10 numéros d'identification appropriés sur 60 molécules A. Coloration de la trajectoire de diffusion de ces 10 molécules. Certes, le mouvement est limité dans la sphère. Je comprends.

Structure et réaction

Enfin, permettez-moi d'expliquer la migration des molécules entre les structures. Premièrement, une espèce moléculaire pour laquelle l'attribut «emplacement» n'est pas spécifié n'appartient à aucune structure. Dans l'exemple précédent, spécifiez «emplacement». Sinon, il sera placé à l'extérieur de la sphère.

Prenons un plan comme exemple et essayons-le. Premièrement, pour faire un plan, nous pouvons le faire comme suit en utilisant trois Real3, l'origine "origine" et les deux vecteurs axes "unité0" et "unité1".

ps = PlanarSurface(origin, unit0, unit1)

En utilisant cela, en supposant un type moléculaire plan A et un type moléculaire normal 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'))

C'est un peu difficile à comprendre, mais les molécules B ne sont disposées que sur le plan. Alors, comment exprimer la réaction dans laquelle ces molécules B sont adsorbées sur le plan M et deviennent des molécules A?

with reaction_rules():
    B + M == A | (1e-3, 1.5)

Comme vous pouvez le voir, lorsque la molécule B entre en collision avec la structure M, elle est adsorbée et devient la molécule A. Au contraire, la molécule A devient le plan d'origine et la molécule B lorsqu'elle est séparée.

À ce stade, vous pouvez calculer avec Simulator comme d'habitude.

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'))

En passant, l'espèce moléculaire de la structure peut être omise en la séparant de la structure, mais elle ne peut pas être omise par liaison.Il est impossible de faire A à partir de B dans la situation où il n'y a pas de M autour duquel la molécule de produit A doit être placée. C'est pourquoi, au contraire, dans ce cas, il n'y a pas de problème car la molécule A peut être placée n'importe où sur le plan M. De plus, dans le cas d'une séparation, la structure est primaire ou non omise. C'est toujours une réaction, mais lors de la liaison, la réaction secondaire devient une réaction primaire et la signification de la constante de vitesse change également.

with reaction_rules():
    B + M > A | 1e-3
    A > B | 3.0  # means the same as A > B + M

Des exercices

1. Exercices sur la méthode spatiale de Gillespie: La méthode spatiale de Gillespie a été utilisée pour calculer l'échelle de temps de diffusion pour un sous-volume de côté $ l $. En revanche, l'échelle de temps nécessaire aux molécules immédiatement après séparation pour se rattacher a également été expliquée. Considérons la dépendance de ces deux échelles de temps sur $ l $ et clarifiez les conditions que $ l $ doit satisfaire, mais on sait que la méthode spatiale de Gillespie devrait satisfaire les conditions suivantes pour $ l $. Oui.

R^2 \ll l^2 \ll 6D\tau_\mathrm{min}

 Où $ R $ est le diamètre de la molécule (plus précisément, le rayon de réaction) et $ \ tau_ \ mathrm {min} $ est la durée de vie la plus courte pour chaque réaction dans la molécule [^ 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. Exercice sur l'arrangement initial des molécules: Modifiez l'exemple ci-dessus pour étudier comment l'arrangement initial des molécules AB, les constantes de vitesse de réaction, les coefficients de diffusion et d'autres paramètres affectent la réaction. Est-il possible de créer un modèle qui reproduit la série chronologique comme l'exemple ci-dessus en utilisant la méthode uniforme de Gillespie?

 1. Exercice sur la diffusion moléculaire: Dans le cas de la diffusion libre en trois dimensions, le déplacement carré moyen est la distance quadratique moyenne à la position de la molécule avec un coefficient de diffusion de $ D $ après $ t $ secondes.

    ```math
\left<L^2\right>=6Dt

Est-ce que c'est effectivement le cas avec Spatiocyte?

1. Exercices sur les constantes de vitesse macroscopiques et microscopiques: les constantes de vitesse macroscopiques et microscopiques ont la même unité de $ \ mathrm {nM} ^ {-1} s ^ {-1} $ et $ \ mathrm {\ mu m} ^ 3s ^ {-1} $. Deux molécules avec un taux de diffusion de 1 $ \ mathrm {\ mu m} ^ 2 / s $ et un diamètre de 5 $ \ mathrm {nm} $ Quelle est la constante du taux de collision $ k_D $ $ \ mathrm {nM} ^ {-1} s ^ {-1} $?

2. Exercice sur la structure: A l'aide d'une structure plane, vérifier la relation entre le coefficient de diffusion $ D $ sur le plan et sa distance quadratique moyenne.