可在www.sciencedirect.com在线获取理论计算机科学电子笔记313（2015）3-21www.elsevier.com/locate/entcsHSIM：一个系统生物学混合随机模拟系统帕特里克·阿马尔1LRI，UMR CNRS 8623 INRIA AMIBUniversit′eParisSud，91405Orsaycedex，法国LoıcPaulev'e2LIX，E´colePolytechnique，91128PalaiseauCedex，FranceETH Zurich，Physikstrasse 3，8092 Zurich，Switzerland摘要本文介绍了一种方法，以获得在一定意义上最好的两个众所周知的模拟方法：化学反应的随机模拟和智能中心（多代理）系统。随机模拟算法及其增强版本适用于模拟充分搅拌的反应物混合物的行为时，空间定位是不相关的;相反，在实体为中心的方法，由于空间定位的原因在计算效率方面，SSA方法受反应数的限制，而不是受反应物的量的限制，而实体中心方法则相反。我们将展示如何从每种方法中受益，并提出一种新的近似随机模拟优化算法。关键词：随机模拟，随机中心系统，混合系统。1引言许多生物对象（遗传调控网络、代谢网络、信号通路等）生物学家和建模者感兴趣的是，不需要考虑所研究过程中涉及的生物分子的空间位置。在有大量分子的情况下，即使是随机方法也不是绝对必要的，而确定性方法，如普通微分方程，给出了有用的结果。但相反，许多生物过程涉及少量的某些分子种类。在一些1电子邮件：pa@lri.fr2电子邮件：loic. ens-cachan.orghttp://dx.doi.org/10.1016/j.entcs.2015.04.0161571-0661/© 2015作者。出版社：Elsevier B.V.这是一篇基于CC BY-NC-ND许可证的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）。4P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3在某些情况下，这些分子或一些更丰富的物种的空间定位是整个系统行为的关键原因在大肠杆菌中，lac操纵子的四聚体阻遏物。大肠杆菌就是一个很好的例子[8]。在每个细菌中只有大约10个拷贝的这种阻遏物，但由于它在DNA上的结合位点附近的局部浓度很高，它的效果在很大程度上是明显的。这种局部高浓度的原因是这种阻遏物的运作机制所固有的，也是它在DNA上结合位点之间的特定距离所固有的这是一个合作的过程，往往会使自己永久化。为了研究这些生物过程，我们需要一个有充分依据的框架，既可以解释随机性和参与者数量少，也可以解释分子的空间定位。最后，它的实现必须在计算时间方面是有效的我们在这里提出的HSIM模拟器，实现了一个混合模拟混合实体为中心的模型与一个新的近似随机模拟Al-出租m（SSA）。HSIM已成功应用于研究涉及空间特征的蛋白质大组装体[1，2]。相关工作。[7]最近介绍了一种混合空间模拟器，涉及精确SSA和单个粒子跟踪。与HSIM相反，他们的方法允许通过切换到单个颗粒或均匀分布的分子来动态改变对物种的处理。2权力中心模型以实体为中心的方法类似于多智能体系统：每个分子都是放置在环境中的一种智能体;计算机根据其特征（化学物种），与其他智能体的相互作用（生化反应）和环境（主要是界定隔室的膜）计算时间和空间的轨迹与多智能体系统的主要区别在于，所有类型的实体都具有相同的行为，当然，这些行为由它们的类型参数化，但基本上是相同的。我们将受益于这种简化，以获得系统的有效这种常见行为有两个组成部分：(i) 扩散：每个分子根据布朗运动扩散到隔室中，其速度取决于温度、介质的粘度和分子拥挤。(ii) 反应：由于它的分散性和与另一个适当物质的分子发生足够能量的碰撞，可以发生反应，导致复合物的形成或解离、消失或一个或两个分子的类型改变P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）35τi/=ji/=j2.1狄氏运动与布朗运动如果我们观察水中的小颗粒，我们会发现它们的运动是不规则的。这种行为，布朗运动，是由于与水分子的碰撞，这使得粒子进行随机行走运动。如果我们观察这样一个粒子的时间间隔为t，我们将看到它的位移为位置向量r。如果我们一遍又一遍地重复这个过程，我们可以计算时间间隔t内的平均位移，我们应该发现= 0。这就意味着，对于随机力，粒子上的水分子布朗运动的特征在于平均位移平方：< |R|2>= 6 Dt其中，D是微分常数。扩散常数D既取决于颗粒的大小和形状，也取决于它所穿过的流体的粘度和温度。因子6与维度有关，6用于3D运动（4用于2D，2用于线性运动）。基本上，微分是由于许多分子的独立布朗运动。2.2随机游走这里有一个简单的模型，它表明3D空间中的随机行走粒子服从与布朗运动相同的平均位移平方定律让则N步后的位移为：r= l−n→1+ l−n→2+. . . +l−n→N其中−→ni是随机方向上的单位向量。因此，N步后粒子所移动的距离的平方为：|R|2=l2（−n→1+−n→2+. . . +−n→N）2|2 = l 2|−n→1|2个以上|−n →2|2+ ..|2+. . . +的|−n→N|2+π−→ni. −n→j=Nl2+l2−→ni. −n→j现在，如果我们重复多次这种随机行走，我们会发现两个不同步骤的平均点产量为零，<−→ni。-n→j>=0，因为步长是随机的方向所以N步后粒子的平均位移平方是< |R|2>= N12. N步后的时间正好是t=Nτ，因此位移平方随时间线性增长：L2< |R|2>=t这意味着这种随机游动是一种布朗运动的模型，其离散常数D=12/ 6τ。6P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3因此，在我们的以实体为中心的模拟器中，离散部分将使用分子的随机行走来进行。将选择参数以说明在真实细胞中体内观察到的蛋白质和小分子的扩散速度2.3反应反应可以分为两大类：单分子和双分子。涉及两种以上反应物的反应可以简化为双分子反应的组合。单分子反应涉及一个可以以多种不同状态存在的分子。蛋白质结构从一种构象到另一种构象的自发波动是单分子反应的一个例子另一个例子是大分子复合物的解离。在以实体为中心的随机模拟器中，这些单分子反应将被建模为具有被触发的概率。为了触发双分子反应，必须发生两个事件：首先，每个反应物分子的副本必须碰撞，其次，这种碰撞必须具有足够的能量。碰撞是由于分子的扩散而产生的;不难看出，碰撞的平均数与反应物浓度的乘积成正比。这就是所谓的群众行动法，由古德伯格和瓦格在1864年提出。这将在模拟器中通过在扩散阶段测试分子的随机移动是否导致与另一个分子的碰撞来实现。如果是，则根据代表其动力学的概率触发由于分子是实体，因此模拟器可以使用这些实体之间的物理链接来考虑组件，从而减缓它们的差异。此外，大的聚集体会使空间变得杂乱，可能导致其他分子的扩散率发生显著变化，从而导致异常的局部浓度（更高或更低）。反过来，这些奇点可以决定整个系统的行为。这些效应中没有一个可以用更全局的方法来考虑，这些方法假设分子是永久均匀的混合物。这一特性带来了计算时间上的代价，即使在没有拓扑效应出现的情况下，全局方法也会更有效。3混合动力系统HSIMHSIM是Hyperstructure SIMulator的缩写，是一个模拟系统，最初开发时考虑了称为Hyperstructures的大型分子组装体的聚集和解离。这些超结构使细菌细胞得以结构化，它们可能在代谢调节、细胞分裂、细胞分化等方面起关键作用。模拟系统必须考虑复合物形成P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）37/解离，它们在细胞中的定位，它们对代谢和信号传导途径以及对调控基因网络的作用。3.1形状、大小和时间尺度HSIM主要用于细菌建模，但也可用于酵母或真菌等小型真核细胞。主隔室是一个圆柱体，两端各有一个半球它可以像许多细菌一样被设置为杆状子隔室可以嵌套到主隔室中。细胞的大小大约是微米。E.大肠杆菌细胞为小杆状，通常约1 μ × 3 μ。酵母细胞是直径约2μ的小球。大多数人类细胞直径约为20μ，细胞核直径约为10μ用来描述分子的基本长度尺度是纳米。许多生物分子是由化学单位（例如核苷酸，氨基酸）组成的，其大小约为1 nm。其他小分子，如脂质和糖，也大致在这个长度范围内。大多数折叠的蛋白质直径只有几纳米在细菌胞质溶胶中，由于分子拥挤，扩散常数比在水中低得多。根据实际细胞中观察到的分子平均速度，调整了HSIM中随机游走的时间步长（100μs）和跳跃距离（10nm）3.2内部HSIM模拟器HSIM是由分子之间的反应规则驱动的随机自动机。从本质上讲，每个分子都由一个记录表示，其中包括其类型，位置，大小以及与某些其他分子的链接列表。HSIM从计算机的角度实时跟踪每个分子。基本原理是时间被切成相同持续时间的连续步骤（迭代），并且在每次迭代中，规则被应用于每个分子。这些规则模拟了真实系统中分子之间8P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3选择来源分子S随机希望在本版单分子硫反应选择随机单位中L方向以S为应用反应根据其概率有什么沿L？是（碰撞）没有S可以在那里方向根据它的扩散概率设T为最近分子是的可能的反应在S和T之间？没有要模拟的模型被指定为一组概率重写规则。这些规则具有以下通用形式：A+B→C+D[p（R）]这条规则意味着当MA和MB（两个物质A和B的分子）碰撞时，它们触发了将MA转化为MC和MB转化为MD（分别是物质C和D）的反应，概率为p（R）。该概率反映了碰撞具有足够的能量来引发反应的事实，并且与反应动力学的宏观解释规则左边的第二个分子同样地，规则右边的第二个分子的省略，说明一个反应只产生一个产物。3.3计算算法在HSIM中，空间是连续的，因为分子坐标是用实数编码的。时间是离散的，它被切割成相等持续时间的片段，在此期间进行一步离散，并引发由此产生的反应。假设只涉及实体，这里是用于计算系统演化的一个步骤的算法的概述Fig. 1. 在每个分子的每个时间片执行算法。在每一步，所有的分子都被认为是一个且只有一次，以随机顺序，以避免统计偏差，因此正确地模拟每个实体的一个过程。分子S的检查包括首先验证它在同一迭代期间没有被预先处理，然后选择随机取向的矢量。的。如果沿着这个矢量没有障碍物（另一个分子或膜），则分子将以与以下相关的概率在其末端移动：它的扩散速度，系统会尝试应用单分子反应。相反，如果与另一个分子T发生碰撞，P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）39在此迭代期间已经检查过，系统将尝试在S和T之间应用双分子反应规则。然后T被标记为“已处理”，并且随机选择下一个分子进行处理，依此类推，直到耗尽所有分子。最后，时间提前了一个时间步，下一次迭代可以开始了（图2）。①的人。如前所述，HSIM是一个混合模拟器：除了由单个实体处理的分子外，它还包括以全局方式处理的其他化学物质这通常是具有高扩散速度的小分子，其以高拷贝数存在，因此在统计学上均匀地分布在隔室内。我们称之为非实体的这些分子类别，将仅由它们在隔室中的类型和拷贝数来表示。使用类似于Gillespie随机方法的全局方法来处理它们将大大提高效率。这个新功能需要在两点上对以前的算法进行轻微修改：(i) 在扩散阶段，如果一个实体可能与非实体物种发生反应，我们必须计算它们在时间片内的平均碰撞次数，然后应用所产生的反应。(ii) 在每次迭代中，在处理分子-实体之前或之后，我们必须处理非实体分子之间的特定相互作用：计算它们在时间片内碰撞的平均次数，然后应用产生的反应。人们可以注意到，使用这种方法，实体之间、实体与非实体之间以及最终非实体之间都可能发生反应。此外，用于描述反应的语法保持相同，无论涉及的分子的类别如何。3.4实体与非实体为了计算在时间片期间实体和非实体之间的平均碰撞次数，计算实体分子所置换的体积，然后，知道非实体分子的浓度，推断它们在该体积内的平均拷贝数是足够的。因为根据定义，非实体物种均匀分布在隔室中，我们所需要做的就是将该特定分子物种的总拷贝数乘以实体取代的体积与隔室总体积的比率。然后，对于这些碰撞中的每一个，根据它们的概率应用相应的反应规则。我们将更准确地看到这是如何实施的在下一节的HSIM3.5一个简单的例子为了演示HSIM提出的三种模拟方法，我们使用了两种酶E1和E2级联的简单例子，其中第一种酶的产物S2是第二种酶的底物。为了使这个例子简单，两种酶的动力学是相同的。为了强调随机效应10P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3在模拟中，代谢物的浓度非常低：对于1微米立方的隔室体积为100个拷贝（0.166μM）。E1+S1→e1s1[0分。1];E2+S2→e2s2[0分。1];e1s1→E1+S1[0分。001];e2s2→E2+S2[0分。001];e1s1→E1+S2[0分。01];e2s2→E2+S3[0分。01];第一种方法是通过指定所有分子都用实体建模来说明的;图2A中的曲线显示了代谢物s1，s2和s3的数量演变的轨迹。A BC D图二. A：使用HSIM实体方法模拟两酶级联模型的轨迹，B：对酶使用HSIM实体，对代谢物使用非实体，C：使用HSIMSSA方法，D：使用SCILABODE分辨率，为了说明第二种方法，酶被建模为实体，而底物被全局处理;图2B中的曲线显示了代谢物s1，s2和s3的数量演变的轨迹。最后，第三种方法是通过指定所有分子都要全局处理来说明的，使用HSIMSSA;相应的曲线如图2C所示。为了了解反应物浓度较大且忽略随机效应时的极限曲线，相应的模型如图2D所示。P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3114近似随机模拟为了 精确模拟 大量不 同物种 的混合物 的进化 ，我们 通常使用 随机模 拟算法（SSA），例如Gillespie的精确直接方法[5]。精确SSA算法的主要缺点之一是当要模拟大量的快速反应时，其计算成本高为了克服这个限制，已经绘制了各种近似SSA，其中大多数是基于τ-跳跃方法[6]。tau-leap方法不是计算下一个发生的反应并相应地更新时间，而是估计在固定时间间隔（τ）内发生的反应数量并立即应用它们这样的算法忽略了应用的反应之间的干扰，因此提供了一个近似的系统的演变。HSIM实现了一个类似于tau-leap方法的近似SSA算法它的主要独创性在于评估在固定时间间隔内发生的反应数量，这需要生成非常少的（均匀的）随机变量。因此，它的总复杂度非常低，其结果是减少了大量反应的时间方差（但保留了它们的平均持续时间）。4.1算法正如前一节所介绍的，HSIM中的每个反应规则（一元或二元）都与在固定时间步长内被触发的概率相关联。在本节的其余部分，我们用p（R）表示反应R的这样一个概率。HSIMSSA的主要仿真循环如下：在每个时间步长τ，（i）计算可能反应的数目，（ii）估计有效反应的数目（根据p（R）），（iii）立即应用这些反应。如果我们只计算有效反应的平均数，结果将接近于常微分方程积分。为了引入一些随机性，反应的（整数）计数被随机舍入：给定一个实数x，将其小数部分与区间[0; 1]（记为U[0;1]）中均匀分布的随机数进行比较，以决定x是舍入到下整数还是上整数。[x]表示x的整数部分。randround（x）：(i) i← [x](ii) 如果U[0;1]≤x−i，则i←i+1(iii) 雷切岛在时间间隔τ期间并且对于每个反应R，HSIM将碰撞（R）评估为可能发生的次数R。在一元反应（1）的情况下，它只是转化物种的种群（[A]是物种A的种群）。在二元反应（2）的情况下，它是可能的碰撞总数，按比例缩放12P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3由细胞体积V和常数α。碰撞（A→C）=Δ[A]（1）碰撞（A+B→C+D）=Δα。[A][B]if（二）V[A]（[A]− 1）/2，否则。然后，通过随机舍入碰撞数量乘以这些碰撞触发反应的概率来完成有效反应R的计数。计数（R）：(i) s←碰撞（R）(ii) n← randround（s. p（R））(iii) if（Ris binary）n←min（n，[RA]，[RB]）(iv) 返回;返回其中[RA]（resp. [B]）代表第一个（分别）的人口。第二）反应R的反应物（见4.3小节）。最后，在每个恒定的时间步长τ，主模拟循环以随机顺序应用每个反应，次数与其有效次数一样多（由count（R）函数返回）。可以通过以下步骤进行总结，其中t是总模拟时间：(i) t←0(ii) 永远重复(a) 对于随机排序的每个反应Rn←count（R）.将反应式R进行n次。(b) t←t+τ总的来说，一个时间步的计算复杂度与反应的数量是线性的，并且包括为每个反应生成一个均匀的随机变量。此外，反应表必须随机排列，这可以在线性时间内完成，R。总的复杂度是O（#R·t），其中#R是反应的数量。在当前的H sim实现中，τ固定为100 μ s（见3.1）; α则被计算为7。4· 10−7，以匹配在具有相同物种浓度的实体中心模型4.2HSIM参数与标准反应速率的关系给定反应R，HSIM假设参数p（R）指定反应在时间步长τ内有效的概率。这个参数应该只取决于所涉及的分子的性质和化学反应的种类。另一方面，标准SSAa laGillespie假设每个反应R由速率（记为r（R））参数化，给出在一个时间单位内可以触发该速率嵌入了细胞的体积P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）313以及一次碰撞后引发反应的概率。通过这种方式，HSIM允许指定独立于细胞体积的参数，从而在体积变化时提供透明的速率缩放p（R）和r（R）之间的关系可以直接从HSIMSSA导出，其中时间间隔由τ时间单位组成1一元反应：r（R）=τ1p（R）; 二元反应：r（R）=ταp（R）。V因为p（R）必须在[0; 1]中，所以该关系对反应速率施加了最大边界，这取决于τ、α和单元体积V。假设这些参数是固定的，这揭示了一个确切的SSA之间的权衡，其成本跟随发生的反应的数量（但与绝对模拟时间无关），而HSIMSSA具有其成本跟随总模拟时间，其可能必须重新缩放以满足上述关系。4.3讨论如4.1小节所述，所提出的算法通过随机变量n = randround（s）来近似反应R的激发次数。其中s=冲突（R）是可能应用的计数。这个随机变量n实际上是s的整数部分之和。p（R），记为i，以及成功概率为s的非整数部分的伯努利随机变量。p（R），注f。因此，n的期望值是i + f = s。p（R），其方差为f−f2。经典的τ-跳跃方法通常用泊松随机变量P（a）来近似反应R的激发次数。r（R）.τ），其中a是所谓的反应R的活度：对于一元反应A→B，a= [A]，对于二元反应A+B→C，a= [A][B]，如果a/=B，否则a= [A]（[A]− 1）/2。从上一小节中建立的p（R）和r（R）之间的关系可以清楚地看出，a。r（R）.τ = s. p（R）。因此，HSIM近似SSA具有与经典τ -跳跃方法相同的期望值，但具有低得多的方差，特别是当i（s的整数部分）时。p（R）是大的。值得注意的是，如上所述，HSIMSSA算法不能经历负种群。一个种群只有在两种情况下才能变成负数(i) 假设有两个反应，Ri和Rj，它们消耗物质A。当应用反应Ri的结果导致消耗足够的反应物A的群体，使得应用反应Ri导致消耗比反应物A的剩余群体更多的反应物A时。(ii) 当[B]/[A]> V/α时，反应A + B → C。p（R），其中[B][A]，最坏的情况是[A] = 1。这可能导致反应的计算应用数量大于其反应物之一的群体在HSIMSSA算法中，第一种情况永远不会出现，因为反应的应用是按顺序进行的，下一次调用计数将使用已经14P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3更新了反应物的数量。这就是为什么为了避免统计偏差，在每个时间步，反应在应用之前随机排序第二种情况在HSIMSSA算法中通过在二元反应的情况 我们可以注意到，只有当一种反应物的数量比另一种反应物的数量大得多时，这种特殊情况才会出现。使用H sim的当前设置（α = 7. 4·10 −7），这导致在生理条件下很少达到的非常高的浓度（几个毫磨牙的数量级）。4.4基本示例我们比较了HSIMSSA算法与两个众所周知的方法：确切的Gillespie SSA方法实现在BIOCHAM和ODE决议的方法实现在SCILAB。对于后者，利用HSIM的输出功能自动获得离散方程组。为了实现这一目标，我们使用了三个基本的例子：(i) 单分子反应将物质A以给定速率转化为物质B的单一(ii) 二元模型，其中双分子反应降解物质A和B的对。(iii) 在种群A和种群B之间存在某种竞争条件的并发模型4.4.1单一模型该模型使用HSIM建模语言显示如下。它简单地说，有两个不同的物种A和B，它们可以在一个隔室中分开。一个微米立方体积。模型的唯一反应是以10−5的概率将A的每个拷贝转换为B的一个拷贝。最初，在隔间内有1000个物种A的拷贝。几何形状= 142：110;// 1 μ ^3隔室。代谢物A、B;A-> B[1e-5];int n（n）;图3A和3B示出了使用HSIMSSA和精确Gillespie方法的模型的一个样本轨迹;可以看出，这些曲线非常相似。图3C中的这些完美的曲线可以通过绘制多次运行的HSIM的平均值来近似（图3D）。误差条显示每个点的标准偏差。P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）315A BC D图三. A：使用HSIMSSA方法的单一模型的模拟轨迹，B：使用BIOCHAMSSA精确方法，C：使用SCILABODE分辨率，D：使用20次HSIM运行获得的平均值和标准差。4.4.2二元模型在这个模型中，有两个不同的物种A和B，它们可以在一个微米立方体积的组件中扩散。该模型的唯一反应是将每对A型和B型分子转化为单个C型分子，概率为0。135.物质C的寿命非常短：如代谢产物C[1]声明所示，该物质的每个拷贝在1个时间步后消失，因此基本上该反应降解了A和B对。最初，在隔室中存在物种A和B中的每一个的1000个几何形状= 142：110;// 1 μ ^3隔室。代谢物A、B、C[1];return（A，B）;A + B-> C[0.135];int n（n）;int n（n）;同样，图4A和4B显示了使用HSIMSSA和精确Gillespie方法的模型的一个样本轨迹;这些曲线也非常相似。图4C中的曲线显示了相同模型的ODE分辨率以及相应的参数值。图4D）显示了每个点的平均值的绘图，16P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3从多次运行的HSIM中获得。A BC D见图4。A：使用HSIMSSA方法的二元模型模拟轨迹，B：使用BioCHAMSSA精确方法，C：使用SCILABODE分辨率，D：使用20次HSIM运行获得的平均值和标准偏差。4.4.3concur模型在这个模型中，物种A在两个未来之间竞争：它可以消失（A→D）或通过与自身反应创造物种B（A+A→B）。同时，物种B也有两种可能的未来：合成两个A的拷贝（B→A+A）或不可逆地转化为C（B→C）。 与前面的例子一样，物质D的寿命很短，所以反应A→D降解A的概率为10−4。 最初在隔间里有10000个物种A的拷贝。代谢物A、B、C、D[1];A->D[1 e-4];A + A-> B[0.2705];B-> A + A[5e-5];B-> C[4e-6];init（10000，A）;该模型的行为非常依赖于反应的动力学在很短的时间内，大约80%的物种AP. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）317A BC D图五. A：使用HSIMSSA方法的并发模型的模拟轨迹，B：使用BioCHAMSSA精确方法，C：使用SCILABODE分辨率，D：使用20次HSIM运行获得的平均值和标准差。转化为B，而其余部分被降解。然后大部分B被转换成C，而剩下的B被转换回A。与前面的例子一样，不同的方法显示了相似的曲线。我们可以注意到，在图5D中，物种A的下降并不像其他模拟中那样快。这是由于用于使误差条可见的量化水平较低。5案例研究5.1蓝藻中的生物钟蓝细菌Synechococcus elongatus是一种革兰氏阴性细菌，利用光合作用获取能量。它的新陈代谢是由一个内部的生物钟通过光敏受体与日光同步调节的。生物钟机制主要由三种蛋白质KaiA、KaiB和KaiC组成。在蓝细菌中，这些蛋白质受到遗传调节，但在体外KaiA、KaiB、KaiC和ATP的混合物中观察到KaiC六聚体的磷酸化水平的自发降解，显示约24小时的时间段[9]。18P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3Kaia磷酸化紧张别构转变别构转变放松去磷酸凯卜Kaia磷酸化T0T1T2T3T4T5T6洗牌R0R1R2R3R4R5R6去磷酸解释这种生物钟运作的模型来自于K. Egushiet al. [4]. KaiC组装成六聚体，可以处于两种状态，紧张或放松（图2）。6A）。A B图第六章 答：KaiC六聚体状态转换和磷酸化水平由KaiA和KaiB调节乙：解释同步振荡的模型一个非常缓慢的变构转换可以将六聚体从一种状态转换到另一种状态。 KaiA蛋白在紧张状态下激活KaiC六聚体的磷酸化。在松弛状态下，KaiC六聚体倾向于去磷酸化。这种去磷酸化被KaiA抵消，KaiA本身被Kaib抑制。KaiA、KaiB和KaiC蛋白的相对浓度导致KaiC六聚体在紧张状态下的磷酸化和在松弛状态下的去磷酸化。在体外实验中，磷酸化和去磷酸化的动力学-KaiC的调节是恒定的并且具有相同的数量级（尽管由于KaiA和KaiB的调节而在体内可变）。所提出的模型显示了每个KaiC六聚体穿过的循环，在此期间其磷酸化水平变化。文献中提出的生物学假设是，单体shu化是所有KaiC六聚体磷酸化水平同步的主要原因，导致系统的稳定周期性行为（图1）。6B）。在该模型中，分子Ri（resp.Ti）代表放松（resp）。紧张）六聚体与I磷酸化。规则（如下所述）包括松弛的单体的舒缩，以及R6与张力单体的舒缩诱导六聚体的松弛。概率已经被固定到Ps=0。5，Pp=10−5，Pr=10−8，Pt=5。10-6i，j，k，l ∈ {0，.，6} s. t. k=i−1，l=j+ 1，Ri+Rji ∈ {0，.，6}，Ti+R6→Rk+Rl→Ri+R6[Ps][Ps]i ∈ {0，.，5}，Ti+R6→Ri+1+R5[Ps]i ∈ {0，.，5}，Ti→Ti+1[Pp]i ∈ {1，.，6}，Ri→Ri−1[Pp]T6→R6[Pr]R6→T6[Pt]P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3192323A B图第七章 A：使用HSIMSSA的KaiC模型的模拟跟踪B：使用BIOCHAMSSA（tau-leap）。该模型的实施导致了非常大量的快速反应，主要是由于Ri单体的收缩。图7A示出了利用HSIM获得的仿真的输出。模拟在1000个单位时间后停止（图中仅显示450个单位），在标准台式计算机上大约花费6秒。为了进行比较，使用近似tau-leap方法（图3）与BioCHAM [37B）需要大约55秒（确切的SSA大约20分钟）。5.2肌动蛋白链聚合为了强调分子的空间定位的重要性，我们提出了一个简化的模型，在一个细胞室中的聚合动力学。这个例子来自真核细胞的肌动蛋白细胞骨架。该隔室最初填充有球状肌动蛋白，球状肌动蛋白是一种可以自组装形成极化二聚体的蛋白质，然后从初始二聚体的每个末端聚合，其动力学取决于聚合过程从初始二聚体的哪个末端开始这个过程是可逆的，聚合物可以收缩，释放球状肌动蛋白单体，经过一定时间后，可以再次聚合。这种机制需要能量，这是由ATP转化为ADP提供的;在这个简化的模型中，我们没有考虑它。肌动蛋白聚合成直线状的纤维，具有缓慢的扩散速度。在这个模型中，我们假设它们根本不扩散。因此，初始二聚体的方向决定了最终薄膜的方向该模型仅用三对规则来定义（“M A M B“意味着两个分子是结合的，而 是直线的）(1) Actg+Actg→MePe[P f]MePe→Actg+Actg[P r]1(2) Actg+Pe→Pe=Actf[P f](3) Actg+Me→Me=Actf[P f]1PeActf→Actg+Pe[P r]MeActf→Actg+Me[P r]第一对规则是二聚体的形成/解离，第二对是从正端聚合/解聚，最后一对是20P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）3负端也是一样。为了解释从正端开始的有利动力学，概率为：Pf/Pr>Pf/Pr。 所有的离散速度2 2 3 3电影的元素（Me，Pe和Actf）被设置为0。在短时间内，游离和结合肌动蛋白的浓度达到平衡并保持稳定。令人惊讶的是，这些元素的几何排列方式：它们似乎或多或少地沿着隔间的大轴线排列（图1）。（八）。图八、 肌动蛋白纤维沿着主隔室的大轴自我排列这种新出现的现象是由于对空头施加的选择性压力。短纤维更倾向于不抵抗长时间的解聚，它们只是消失。相反，在同样的条件下，长链可能会变短，但当它们保持不变时，它们的取向也保持不变。因此，剩下的纤维沿着隔间的大轴，在有足够空间生长的地方当然，如果没有一个考虑到空间定位的系统， 和分子集合体的几何形状，我们就不会意识到这种现象。6结论HSIM实现了一个混合仿真系统，它集成了以实体为中心的系统和近似随机仿真算法（SSA）。我们已经展示了以实体为中心的仿真系统在表现力方面的优势，特别是考虑到空间定位。我们还提出了一个新的和非常有效的近似SSA计算系统的反应规则的演变。我们通过两个非常不同的生物模型的研究来说明这种方法。第一个不需要空间定位，但涉及大量反应，与标准近似SSA相比，HSIM所需的计算时间减少了近10倍。第二个，相反，空间定位和考虑HSIM分子聚集体的能力，已被证明是发现一个涌现现象的决定性因素。P. 阿马尔湖Paulevé/Electronic Notes in Theoretical Computer Science 313（2015）321H SIM可以在http://www.lri.fr/ ~ pa/Hsim下载。引用[1] Amar P.，Bernot G.，诺里斯五世，HSIM：一个模拟程序来研究蛋白质的大组件，生物物理学和化学杂志4（2004），79-84。[2] Amar P.，Legent G.，Thellier M.，Ripoll C.，伯诺特·G Nystrom T.，小赛尔，诺里斯五世，随机自动机显示酶组装如何有助于代谢效率，BMC Systems Biology2（2008），27。[3] Calzone L.，Fages F.，Soliman S.，BIOCHAM：模拟生物系统和形式化实验知识的环境，生物信息学，22（2006），1805-1807[4] 江口，尤达·M寺田T.P.，和Sasai M.，KaiC体外磷酸化的稳健昼夜节律振荡机制，生物物理学杂志95（4）（2008），1773-1784。[5] Gillespie ，Daniel T.，耦合化学反应的精确随机模拟，物理化学杂志81（25）（1977），2340-2361。[6] Gillespie，Daniel T.，化学反应系统的近似加速随机模拟，物理化学杂志115（2001），1716-1733。[7] Klann M. ， Ganguly A. ， 和 Koeppl H. ， 混 合 空 间 Gillespie 和 粒 子 跟 踪 模 拟 ， BMC Bioinformatics28（2012），18，pp.i549-i555。[8] Müller-Hil，Benno，作为辅助lac算子的距离、相位和质量的函数的lac算子的表示。J. Mol. 257（1996），21-29。[9] Nakajima M.，今井，伊藤H西胁T.，村山，岩崎H.，小山T.，和近藤T.，Reconstitution of CircadianOscillation of Cyanobacterial KaiC Phosphorylation in Vitro，Science308（2005），414-415.
我的内容管理
展开
