如果我使用linked question 中的可运行代码并替换您对initial_guess 的定义：
mean_gauss_x = sum(x * data.reshape(201,201)) / sum(data.reshape(201,201))
sigma_gauss_x = np.sqrt(sum(data.reshape(201,201) * (x - mean_gauss_x)**2) / sum(data.reshape(201,201)))
mean_gauss_y = sum(y * data.reshape(201,201)) / sum(data.reshape(201,201))
sigma_gauss_y = np.sqrt(sum(data.reshape(201,201) * (y - mean_gauss_y)**2) / sum(data.reshape(201,201)))
initial_guess = (np.max(data), mean_gauss_x, mean_gauss_y, sigma_gauss_x, sigma_gauss_y,0,10)
然后
print(inital_guess)
产量
(13.0, array([...]), array([...]), array([...]), array([...]), 0, 10)
注意initial_guess 中的一些值是数组。 optimize.curve_fit 函数期望 initial_guess 是一个标量元组。这就是问题的根源。
错误信息
ValueError: setting an array element with a sequence
通常在需要标量值时提供类数组时出现。这暗示问题的根源可能与维数错误的数组有关。例如，如果您将一维数组传递给需要标量的函数，则可能会出现这种情况。
让我们看看这段取自linked question的代码：
x = np.linspace(0, 200, 201)
y = np.linspace(0, 200, 201)
X, Y = np.meshgrid(x, y)
x 和 y 是一维数组，而 X 和 Y 是二维数组。 （我已将所有二维数组大写以帮助将它们与一维数组区分开来）。
现在请注意，Python sum 和 NumPy 的 sum 方法在应用于二维数组时表现不同：
In [146]: sum(X)
Out[146]:
array([ 0., 201., 402., 603., 804., 1005., 1206., 1407.,
1608., 1809., 2010., 2211., 2412., 2613., 2814., 3015.,
...
38592., 38793., 38994., 39195., 39396., 39597., 39798., 39999.,
40200.])
In [147]: X.sum()
Out[147]: 4040100.0
Python sum 函数等价于
total = 0
for item in X:
total += item
由于X 是一个二维数组，循环for item in X 正在遍历X 的行。因此，每个item 都是一个表示一行X 的一维数组。因此，total 最终成为一维数组。
相比之下，X.sum() 将X 中的所有元素相加并返回一个标量。
因为initial_guess 应该是一个标量元组，
在您使用 sum 的任何地方，您都应该使用 NumPy sum 方法。例如，替换
mean_gauss_x = sum(x * data) / sum(data)
与
mean_gauss_x = (X * DATA).sum() / (DATA.sum())
import numpy as np
import scipy.optimize as optimize
import matplotlib.pyplot as plt
# define model function and pass independant variables x and y as a list
def twoD_Gaussian(data, amplitude, xo, yo, sigma_x, sigma_y, theta, offset):
X, Y = data
xo = float(xo)
yo = float(yo)
a = (np.cos(theta) ** 2) / (2 * sigma_x ** 2) + (np.sin(theta) ** 2) / (
2 * sigma_y ** 2
)
b = -(np.sin(2 * theta)) / (4 * sigma_x ** 2) + (np.sin(2 * theta)) / (
4 * sigma_y ** 2
)
c = (np.sin(theta) ** 2) / (2 * sigma_x ** 2) + (np.cos(theta) ** 2) / (
2 * sigma_y ** 2
)
g = offset + amplitude * np.exp(
-(a * ((X - xo) ** 2) + 2 * b * (X - xo) * (Y - yo) + c * ((Y - yo) ** 2))
)
return g.ravel()
# Create x and y indices
x = np.linspace(0, 200, 201)
y = np.linspace(0, 200, 201)
X, Y = np.meshgrid(x, y)
# create data
data = twoD_Gaussian((X, Y), 3, 100, 100, 20, 40, 0, 10)
data_noisy = data + 0.2 * np.random.normal(size=data.shape)
DATA = data.reshape(201, 201)
# add some noise to the data and try to fit the data generated beforehand
mean_gauss_x = (X * DATA).sum() / (DATA.sum())
sigma_gauss_x = np.sqrt((DATA * (X - mean_gauss_x) ** 2).sum() / (DATA.sum()))
mean_gauss_y = (Y * DATA).sum() / (DATA.sum())
sigma_gauss_y = np.sqrt((DATA * (Y - mean_gauss_y) ** 2).sum() / (DATA.sum()))
initial_guess = (
np.max(data),
mean_gauss_x,
mean_gauss_y,
sigma_gauss_x,
sigma_gauss_y,
0,
10,
)
print(initial_guess)
# (13.0, 100.00000000000001, 100.00000000000001, 57.106515650488404, 57.43620227324201, 0, 10)
# initial_guess = (3,100,100,20,40,0,10)
popt, pcov = optimize.curve_fit(twoD_Gaussian, (X, Y), data_noisy, p0=initial_guess)
data_fitted = twoD_Gaussian((X, Y), *popt)
fig, ax = plt.subplots(1, 1)
ax.imshow(
data_noisy.reshape(201, 201),
cmap=plt.cm.jet,
origin="bottom",
extent=(X.min(), X.max(), Y.min(), Y.max()),
)
ax.contour(X, Y, data_fitted.reshape(201, 201), 8, colors="w")
plt.show()