【发布时间】：2010-06-16 14:14:52
【问题描述】：
*于 2010 年 6 月 17 日编辑
我正在尝试了解如何改进我的代码（使其更具 Python 风格）。此外，我有兴趣编写更直观的“条件”来描述生物化学中常见的场景。下面程序中的条件标准我在答案 #2 中已经解释过，但我对代码不满意——它工作得很好，但并不明显，而且对于更复杂的条件场景也不容易实现。欢迎提出想法。欢迎评论/批评。第一次发帖经验@stackoverflow-如果需要，请评论礼仪。
代码生成一个值列表，作为以下练习的解决方案：
“在您选择的编程语言中，实施 Gillespie 的第一反应算法来研究反应 A--->B 的时间行为，其中从 A 到 B 的转变只有在另一个化合物 C 是存在，其中 C 与 D 动态相互转换，如下面的 Petri 网中所示。假设在反应开始时有 100 个 A 分子，1 个 C 分子，并且没有 B 或 D 存在。将 kAB 设置为 0.1 s -1 和 kCD 和 kDC 都为 1.0 s-1。模拟系统在 100 s 内的行为。"
def sim():
# Set the rate constants for all transitions
kAB = 0.1
kCD = 1.0
kDC = 1.0
# Set up the initial state
A = 100
B = 0
C = 1
D = 0
# Set the start and end times
t = 0.0
tEnd = 100.0
print "Time\t", "Transition\t", "A\t", "B\t", "C\t", "D"
# Compute the first interval
transition, interval = transitionData(A, B, C, D, kAB, kCD, kDC)
# Loop until the end time is exceded or no transition can fire any more
while t <= tEnd and transition >= 0:
print t, '\t', transition, '\t', A, '\t', B, '\t', C, '\t', D
t += interval
if transition == 0:
A -= 1
B += 1
if transition == 1:
C -= 1
D += 1
if transition == 2:
C += 1
D -= 1
transition, interval = transitionData(A, B, C, D, kAB, kCD, kDC)
def transitionData(A, B, C, D, kAB, kCD, kDC):
""" Returns nTransition, the number of the firing transition (0: A->B,
1: C->D, 2: D->C), and interval, the interval between the time of
the previous transition and that of the current one. """
RAB = kAB * A * C
RCD = kCD * C
RDC = kDC * D
dt = [-1.0, -1.0, -1.0]
if RAB > 0.0:
dt[0] = -math.log(1.0 - random.random())/RAB
if RCD > 0.0:
dt[1] = -math.log(1.0 - random.random())/RCD
if RDC > 0.0:
dt[2] = -math.log(1.0 - random.random())/RDC
interval = 1e36
transition = -1
for n in range(len(dt)):
if dt[n] > 0.0 and dt[n] < interval:
interval = dt[n]
transition = n
return transition, interval
if __name__ == '__main__':
sim()
【问题讨论】：
-
很难从令人窒息的意识流中解析出你的实际问题。您能否总结一下您的实际问题，以便我们实际回答？
-
这似乎是蒙特卡罗位不正确。如果
dt[2]小于 1e36，则转换将始终等于 2。我承认我不完全理解这种模拟，但那部分对我来说似乎很可疑。似乎dts 应该以某种方式用作权重，以随机确定哪个进程将进行。有随机性，但它似乎是错误的方式。 -
谢谢大家。我还没有看过贾斯汀的蒙特卡罗方面，但我会尽快发表评论。 S. Lott-你是对的，我不喜欢那样做。我会稍微编辑一下我的帖子——哈哈，我刚刚发现了“将返回值传递给调用者”的含义——这是我不熟悉的部分。到目前为止，我已经能够通过使用数组定义函数来摆脱一些简单的模拟游戏sigh。但是，是的-我的帖子的要点是-WTF？我认为我在正确的轨道上。今天晚些时候我会重新发布。
标签： python physics stochastic