scipy.optimize.最小二乘法从Python中使用C++程序时不能给出正确的结果

Question

<>我试图使用^ {CD1>}函数，以最小化C++代码生成的模型。在

下面是一个最小的例子来演示这个问题（taken from the Scipy docs）。在

我使用两个函数，一个是纯Python，另一个调用C++代码生成模型。两个函数在测试时返回相同的结果，但在scipy.optimize.least_squares中使用时则不返回相同的结果。在

import numpy as np
from scipy.optimize import least_squares
import csv
import os
import time

def gen_data(t, a, b, c, noise=0, n_outliers=0, random_state=0):
    y = a + b * np.exp(t * c)
    rnd = np.random.RandomState(random_state)
    error = noise * rnd.randn(t.size)
    outliers = rnd.randint(0, t.size, n_outliers)
    error[outliers] *= 10
    return y + error

a = 0.5
b = 2.0
c = -1
t_min = 0
t_max = 10
n_points = 15
t_train = np.linspace(t_min, t_max, n_points)
y_train = gen_data(t_train, a, b, c, noise=0.1, n_outliers=3)

def fun_cpp(x, t, y):
    a = os.popen("./a.out "+str(x[0])+" "+str(x[1])+" "+str(x[2])).readlines()
    model = [float(i) for i in a[0].split(' ')[:-1]]

    return model - y

def fun_python(x, t, y):
    return (x[0] + x[1] * np.exp(x[2] * t)) - y

x0 = np.array([1.0, 1.0, 0.0])

res_lsq = least_squares(fun_python, x0, args=(t_train, y_train), verbose = 2)
res_lsq2 = least_squares(fun_cpp, x0, args=(t_train, y_train), verbose = 2)

t_test = np.linspace(t_min, t_max, n_points * 10)
y_true = gen_data(t_test, a, b, c)
y_lsq = gen_data(t_test, *res_lsq.x)
y_lsq2 = gen_data(t_test, *res_lsq2.x)

import matplotlib.pyplot as plt
plt.ion() 
plt.plot(t_train, y_train, 'o')
plt.plot(t_test, y_true, 'k', linewidth=2, label='true')
plt.plot(t_test, y_lsq, label='pure Python')
plt.plot(t_test, y_lsq2, label='C++/Python')

plt.legend()
plt.show()

这里是C++代码（运行：^ {CD3>}）：

^{pr2}$

此代码生成下图：

有什么建议吗？在

谢谢

Answer 1

我想我的问题已经解决了。在

{{cd2>不需要调用^ cd2>的外部函数。在

为了解决我的问题，可以使用黑盒优化算法，例如：https://github.com/paulknysh/blackbox

上面的示例代码可以修改为：

import numpy as np
from scipy.optimize import least_squares
import csv
import os
import time

def gen_data(t, a, b, c, noise=0, n_outliers=0, random_state=0):
    y = a + b * np.exp(t * c)
    rnd = np.random.RandomState(random_state)
    error = noise * rnd.randn(t.size)
    outliers = rnd.randint(0, t.size, n_outliers)
    error[outliers] *= 10
    return y + error

a = 0.5
b = 2.0
c = -1
t_min = 0
t_max = 10
n_points = 15
t_train = np.linspace(t_min, t_max, n_points)
y_train = gen_data(t_train, a, b, c, noise=0.1, n_outliers=3)

def fun_cpp(x):
    global y_train
    a = os.popen("./a.out "+str(x[0])+" "+str(x[1])+" "+str(x[2])).readlines()
    model = [float(i) for i in a[0].split(' ')[:-1]]
    chisquared = np.sum([(y1-y2)**2. for y1,y2 in zip(y_train, model)])
    return chisquared

def fun_python(x, t, y):
    return (x[0] + x[1] * np.exp(x[2] * t)) - y

x0 = np.array([1.0, 1.0, 0.0])

res_lsq = least_squares(fun_python, x0, args=(t_train, y_train))

import blackbox as bb

bb.search(f=fun_cpp,  # given function
              box=[[0., 1.], [0., 3], [-2., 1.]],  # range of values for each parameter
              n=40,  # number of function calls on initial stage (global search)
              m=40,  # number of function calls on subsequent stage (local search)
              batch=4,  # number of calls that will be evaluated in parallel
              resfile='output.csv')  # text file where results will be saved

rowSelect = []
with open('output.csv') as csvfile:
    readCSV = csv.reader(csvfile, delimiter=',')
    for row in readCSV:
        rowSelect.append(row)

bestFitParam = [float(x) for x in rowSelect[1][0:-1]]

t_test = np.linspace(t_min, t_max, n_points * 10)
y_true = gen_data(t_test, a, b, c)
y_lsq = gen_data(t_test, *res_lsq.x)
y_lsq2 = gen_data(t_test, *bestFitParam)

import matplotlib.pyplot as plt
plt.ion() 
plt.plot(t_train, y_train, 'o')
plt.plot(t_test, y_true, 'k', linewidth=2, label='true')
plt.plot(t_test, y_lsq, label='pure Python')
plt.plot(t_test, y_lsq2, label='Black-box C++/Python')

plt.legend()
plt.show()

结果如下：

黑盒优化可以获得更好的结果，但计算时间较长。在