我正在尝试改编来自matlab的代码,并想知道如何定义T(用于时间)和X(用于反应选择),以便在我的Gillespie算法运行时对其进行更新
k1 = 0.1
k2 = 1.0
time_limit = 1800
max_length = 1000
count = 0
T[0] = 0
X[0,:] = [100, 0, 0]
N = np.array([[-1, 1, 0],[0, -1, 1]])
X[0,:] = [100, 0, 0]
while T(count,1) < time_limit or count <= max_length:
a = [k1*X[count,0], k2*X[count,1]]
a0 = sum(a)
r = np.random.rand(2)
#step length
tau = -log(r(1))/a0
#reaction choice
cumsuma = np.cumsum(a)
cumsuma = np.cumsum(a)
for i in cumsuma:
if cumsuma[i] >= r(2) * a0:
mu = i
break
T(count+1) = T(count)+tau;
X(count+1,:) = X(count,:) + stoich_matrix(mu,:)
count = count+1
在MATLAB / Octave中,您可以“开始”一个带有索引表达式的变量,并通过索引另一个值使其变大:
>> T(1) = 1
T = 1
>> T(2) = 1
T =
1 1
[我不认为这是我几年前使用MATLAB时的好习惯,但显然允许初学者/懒惰的程序员让步。
您无法在Python / numpy中做到这一点。>>
In [217]: T[0] = 1 --------------------------------------------------------------------------- NameError Traceback (most recent call last) <ipython-input-217-319f5945c3d7> in <module> ----> 1 T[0] = 1 NameError: name 'T' is not defined
您首先初始化了
T
数组:
In [218]: T = np.zeros(10, int) In [219]: T Out[219]: array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0]) In [220]: T[0] = 1 In [221]: T[2] = 3 In [222]: T Out[222]: array([1, 0, 3, 0, 0, 0, 0, 0, 0, 0])
而且您不能将数组用作函数:
In [223]: T(1) --------------------------------------------------------------------------- TypeError Traceback (most recent call last) <ipython-input-223-04be72b98e92> in <module> ----> 1 T(1) TypeError: 'numpy.ndarray' object is not callable
您正确创建了
N
数组:
N = np.array([[-1, 1, 0],[0, -1, 1]])
表达方式:
T(count+1) = T(count)+tau; X(count+1,:) = X(count,:) + stoich_matrix(mu,:)
在Python中是错误的。对于MATLAB,它们看起来还不错,但是
()
用于索引和函数调用。在Python中,foo()
是函数调用,foo[...]
是索引操作。(mu,:)
对于函数调用不利,[mu, :]
对于索引操作是可以的。
如果没有对两种语言的充分了解,很难将其从MATLAB转换为Python。但是,对写作目的有足够的了解是最重要的。