<p>显然要做的是从<code>data</code>中删除nan。然而,这样做还需要移除2D <code>X</code>,<code>Y</code>位置数组中的相应位置:</p>
<pre><code>X, Y = np.indices(data.shape)
mask = ~np.isnan(data)
x = X[mask]
y = Y[mask]
data = data[mask]
</code></pre>
<p>现在您可以使用<code>optimize.leastsq</code>(或更新、更简单的<code>optimize.curve_fit</code>)
要使数据适合模型函数:</p>
^{pr2}$
<hr/>
<p>例如,如果我们用nan生成一些随机的<code>data</code></p>
<pre><code>data = make_data(shape)
</code></pre>
<p>所以</p>
<pre><code>import matplotlib.pyplot as plt
plt.imshow(data)
plt.show()
</code></pre>
<p>看起来像</p>
<p><img src="https://i.stack.imgur.com/gTvxb.png" alt="enter image description here"/></p>
<p>用白点表示有NaN值的地方,然后</p>
<pre><code>import numpy as np
from scipy import optimize
np.set_printoptions(precision=4)
def gaussian(p, x, y):
height, center_x, center_y, width_x, width_y = p
return height*np.exp(-(((center_x-x)/width_x)**2+((center_y-y)/width_y)**2)/2)
def moments(data):
total = np.nansum(data)
X, Y = np.indices(data.shape)
center_x = np.nansum(X*data)/total
center_y = np.nansum(Y*data)/total
row = data[int(center_x), :]
col = data[:, int(center_y)]
width_x = np.nansum(np.sqrt(abs((np.arange(col.size)-center_y)**2*col))
/np.nansum(col))
width_y = np.nansum(np.sqrt(abs((np.arange(row.size)-center_x)**2*row))
/np.nansum(row))
height = np.nanmax(data)
return height, center_x, center_y, width_x, width_y
def errorfunction(p, x, y, data):
return gaussian(p, x, y) - data
def fitgaussian(data):
params = moments(data)
X, Y = np.indices(data.shape)
mask = ~np.isnan(data)
x = X[mask]
y = Y[mask]
data = data[mask]
p, success = optimize.leastsq(errorfunction, params, args=(x, y, data))
return p
def make_data(shape):
h, w = shape
p = 50, h/2.0, w/2.0, h/3.0, w/5.0
print('Actual parameters: {}'.format(np.array(p)))
X, Y = np.indices(shape)
data = gaussian(p, X, Y) + np.random.random(shape)
mask = np.random.random(shape) < 0.3
data[mask] = np.nan
return data
shape = 100, 200
data = make_data(shape)
X, Y = np.indices(shape)
parameters = fitgaussian(data)
print('Fitted parameters: {}'.format(parameters))
fit = gaussian(parameters, X, Y)
</code></pre>
<p>收益率</p>
<pre><code>Actual parameters: [ 50. 50. 100. 33.3333 40. ]
Fitted parameters: [ 50.2908 49.9992 99.9927 33.7039 40.6149]
</code></pre>