Python的双变量搜索算法

Question

我有一个算法，其中：

• 采用两个输入（X，Y），其中，x和y是在两个独立的泊松分布（http://en.wikipedia.org/wiki/Poisson_distribution#Definition）
所述的lambda'变量
• 计算泊松矢量用于每个x和y的，
• 计算两个泊松的矩阵乘积矢量
• 返回[总和（upper_quadrant），总和（lower_quadrant），总和（对角线）]

所以：

import math, random 

def poisson(m, n): 
    p=math.exp(-m) 
    r=[p] 
    for i in range(1, n): 
     p*=m/float(i) 
     r.append(p) 
    return r 

def simulate(mx, my, n): 
    r=[0.0 for i in range(3)] 
    px, py = (poisson(mx, n), 
       poisson(my, n)) 
    for i in range(n): 
     for j in range(n): 
      if i > j: 
       k=0 
      elif i < j: 
       k=1 
      else: 
       k=2 
      r[k]+=px[i]*py[j] 
    return r 

现在，我需要解决的X和Y，给出一组特定的输出。我已经砍死在一起下列求解器功能：

def solve(p, n, generations, decay, tolerance): 
    def rms_error(X, Y): 
     return (sum([(x-y)**2 
        for x, y in zip(X, Y)])/float(len(X)))**0.5 
    def calc_error(x, y, n, target): 
     guess=simulate(x, y, n) 
     return rms_error(target, guess)  
    q=[0, 0] 
    err=calc_error(math.exp(q[0]), math.exp(q[1]), n, p) 
    bestq, besterr = q, err 
    for i in range(generations): 
     if besterr < tolerance: 
      break 
     q=list(bestq) 
     if random.random() < 0.5: 
      j=0 
     else: 
      j=1 
     fac=((generations-i+1)/float(generations))**decay 
     q[j]+=random.gauss(0, 1)*fac 
     err=calc_error(math.exp(q[0]), math.exp(q[1]), n, p) 
     if err < besterr: 
      bestq, besterr = q, err 
      # print (i, bestq, besterr) 
    q, err = [math.exp(q) for q in bestq], besterr 
    return (i, q, err) 

其作品，但似乎需要比较大量的尝试返回一个不是非常优化的响应：

if __name__=="__main__": 
    p, n = [0.5, 0.2, 0.3], 10 
    q, err = solve_match_odds(p, n, 
           generations=1000, 
           decay=2, 
           tolerance=1e-5) 
    print q 
    print simulate_match_odds(q[0], q[1], n) 
    print (i, err) 

和：

justin@justin-ThinkPad-X220:~/work/$ python solve.py 
[0.5, 0.2, 0.3] 
[0.5000246335218251, 0.20006624338256798, 0.29990837191131686] 
(999, 6.680993630511076e-05) 
justin@justin-ThinkPad-X220:~/work/$ 

我不是CS专业，我觉得我错过了这里所有的搜索文献。有人可以提出一个更好的方法来搜索像这样的二维空间中的变量吗？

谢谢。

Answer 1

您可以使用scipy.optimize库来解决这个问题：

import math, random 

def poisson(m, n): 
    p=math.exp(-m) 
    r=[p] 
    for i in range(1, n): 
     p*=m/float(i) 
     r.append(p) 
    return r 

def simulate(mx, my, n): 
    r=[0.0 for i in range(3)] 
    px, py = (poisson(mx, n), 
       poisson(my, n)) 
    for i in range(n): 
     for j in range(n): 
      if i > j: 
       k=0 
      elif i < j: 
       k=1 
      else: 
       k=2 
      r[k]+=px[i]*py[j] 
    return r 

from scipy import optimize 

Target, N = [0.5, 0.2, 0.3], 10 

def error(p, target, n): 
    r = simulate(p[0], p[1], n) 
    return np.sum(np.subtract(target, r)**2) 

r = optimize.fmin(error, (0, 0), args=(Target, N)) 
print simulate(r[0], r[1], n) 

输出：

Optimization terminated successfully. 
     Current function value: 0.000000 
     Iterations: 67 
     Function evaluations: 125 
[0.49999812501285623, 0.20000418001288445, 0.2999969464616799]