PYTHON：Scipy的曲线适合NxM阵列？

问题描述

我通常使用Scipy.Optimize.curveFit来使定制函数适合数据。 这种情况下的数据始终是一维数组。

二维数组有类似的函数吗？

举个例子，我有一个10x10的数值数组。然后我有一个函数，它执行一些操作并创建一个10x10数值数组，我希望对该函数进行拟合，以便得到的10x10数组最适合输入数组。

也许举个例子更好：)

data = pyfits.getdata('data.fits')  #fits is an image format, this gives me a NxM numpy array
mod1 = pyfits.getdata('mod1.fits')
mod2 = pyfits.getdata('mod2.fits')    
mod3 = pyfits.getdata('mod3.fits')

mod1_1D = numpy.ravel(mod1)
mod2_1D = numpy.ravel(mod2)    
mod3_1D = numpy.ravel(mod3)
def dostuff(a,b):    #originaly this is a function for 2D arrays
    newdata = (mod1_1D*12)+(mod2_1D)**a - mod3_1D/b
    return newdata

现在应该安装a和b，以便新数据尽可能接近数据。

到目前为止我得到的：

data1D = numpy.ravel(data)
data_X = numpy.arange(data1D.size)
fit = curve_fit(dostuff,data_X,data1D)

但印刷品只适合我

(array([ 1.]), inf)

我的阵列中确实有一些nan，也许这是一个问题？

---

解决方案

目标是将2D函数表示为一维函数：g(x, y, ...) --> f(xy, ...)

将坐标对(x, y)转换为单个数字xy乍看起来可能很棘手。但实际上这很简单。只需枚举所有数据点，您就有一个唯一定义每个坐标对的数字。拟合后的函数只需重建原始坐标，进行计算并返回结果。

适合20x10图像中的二维线性渐变的示例：

import scipy as sp
import numpy as np
import matplotlib.pyplot as plt

n, m = 10, 20

# noisy example data
x = np.arange(m).reshape(1, m)
y = np.arange(n).reshape(n, 1)
z = x + y * 2 + np.random.randn(n, m) * 3

def f(xy, a, b):
    i = xy // m  # reconstruct y coordinates
    j = xy % m  # reconstruct x coordinates
    out = i * a + j * b
    return out

xy = np.arange(z.size)  # 0 is the top left pixel and 199 is the top right pixel
res = sp.optimize.curve_fit(f, xy, np.ravel(z))

z_est = f(xy, *res[0])
z_est2d = z_est.reshape(n, m)


plt.subplot(2, 1, 1)
plt.plot(np.ravel(z), label='original')
plt.plot(z_est, label='fitted')
plt.legend()

plt.subplot(2, 2, 3)
plt.imshow(z)
plt.xlabel('original')

plt.subplot(2, 2, 4)
plt.imshow(z_est2d)
plt.xlabel('fitted')