代码之家  ›  专栏  ›  技术社区  ›  sun0727

使用scipy得到错误的结果。优化fmin\U cobyla公司

  •  0
  • sun0727  · 技术社区  · 7 年前

    我对scipy非常陌生,现在我很难使用函数 in scipy.optimize 通过做一些小实验。

    我试图通过找到误差值最低的参数来拟合sin函数。

    使用的功能是 fmin_cobyla

    代码如下:

    import matplotlib.pyplot as plt
    from scipy.optimize import fmin_cobyla
    from scipy.optimize import fmin_slsqp
    from scipy.optimize import leastsq
    import numpy as np
    from sympy import *
    
    noise = np.random.randn(100)
    
    def func_model(x, para):
        ''' Model: y = a*sin(2*k*pi*x+theta)'''
        a, k, theta = para
        return a*np.sin(2*k*np.pi*x+theta)
    
    def func_noise(x, para):
        a, k, theta = para
        return a*np.sin(2*k*np.pi*x+theta) + noise
    
    def func_error(para_guess):
        '''error_func'''
        x_seq = np.linspace(-2*np.pi, 0, 100)
        para_fact = [10, 0.34, np.pi/6]
        data = func_noise(x_seq, para_fact)
        error_value  = data - func_model(x_seq, para_guess)
        return error_value
    
    # 1<a<15  0<k<1  0<theta<pi/2
    constraints = [lambda x: 15 - x[0], lambda x: x[0]- 1, \
                   lambda x: 1 - x[1],  lambda x: x[1], \
                   lambda x: np.pi/2 - x[2], lambda x: x[2]]
    
    para_guess_init = np.array([7, 0.2, 0])
    
    solution = fmin_cobyla(func_error, para_guess_init, constraints)
    print(solution)   # supposed to be like [10, 0.34, np.pi/6]
    
    xx = np.linspace(-2*np.pi, 0, 100)
    plt.plot(xx, func_model(xx, [10, 0.34, np.pi/6]), label="raw")
    plt.plot(xx, func_noise(xx, [10, 0.34, np.pi/6]), label="with noise")
    plt.plot(xx, func_model(xx, solution), label="fitted")
    plt.legend()
    plt.show()
    

    跑步后我得到了结果

    解决方案=[1.6655938 0.59868667 0.0731335]

    enter image description here

    这肯定不是正确答案

    谁能帮帮我吗。提前谢谢。。

    1 回复  |  直到 7 年前
        1
  •  3
  •   xnx    7 年前

    这里有两件事看起来显然是错误的:首先,每次调用目标函数时,你都在改变噪声,所以你的优化是试图击中一个移动的目标。调用前设置模拟数据 fmin_cobyla :

    the_noise = np.random.randn(100)
    data = func_noise(x_seq, para_fact)
    

    还有,你的 func_error 应返回模型与每个点的数据之间的差值,而不是平方和差值:

    def func_error(para_guess):
        error_value = data - func_model(x_seq, para_guess)
        return error_value
    

    你还是会发现 fmin\U cobyla公司 努力寻找约束最小值。。。一些预处理可以更好地估计相位或频率的初始猜测,这可能会对您有所帮助。