【问题标题】:Python curve_fit with multiple independent variables具有多个自变量的 Python curve_fit
【发布时间】:2015-04-06 23:57:33
【问题描述】:

Python 的 curve_fit 计算具有单个自变量的函数的最佳拟合参数,但有没有办法使用 curve_fit 或其他方法来拟合具有多个自变量的函数?例如:

def func(x, y, a, b, c):
    return log(a) + b*log(x) + c*log(y)

其中 x 和 y 是自变量,我们希望拟合 a、b 和 c。

【问题讨论】:

    标签: python scipy curve-fitting


    【解决方案1】:

    您可以为自变量传递一个多维数组curve_fit,但是您的func 必须接受同样的东西。例如,调用此数组X 并将其解包为xy 以便清楚起见:

    import numpy as np
    from scipy.optimize import curve_fit
    
    def func(X, a, b, c):
        x,y = X
        return np.log(a) + b*np.log(x) + c*np.log(y)
    
    # some artificially noisy data to fit
    x = np.linspace(0.1,1.1,101)
    y = np.linspace(1.,2., 101)
    a, b, c = 10., 4., 6.
    z = func((x,y), a, b, c) * 1 + np.random.random(101) / 100
    
    # initial guesses for a,b,c:
    p0 = 8., 2., 7.
    print(curve_fit(func, (x,y), z, p0))
    

    适合:

    (array([ 9.99933937,  3.99710083,  6.00875164]), array([[  1.75295644e-03,   9.34724308e-05,  -2.90150983e-04],
       [  9.34724308e-05,   5.09079478e-06,  -1.53939905e-05],
       [ -2.90150983e-04,  -1.53939905e-05,   4.84935731e-05]]))
    

    【讨论】:

    • 如果 x 和 y 的大小不同,是否可以修改解决方案以使用曲线拟合。例如,x = linspace(0.1,1.1,101) 和 y = np.array([1.0,2.0])?
    • 我不确定我是否关注你:func 表示一个二元函数(取两个自变量),因此对于拟合,它应该定义为给出结果 f(x_i,y_i)任何提供的输入值 x_i 和 y_i。如果 xy 的大小不同,那么您正在尝试评估它,例如在一些xy undefined 这肯定不能完成。
    • 我只想分享与@ScottG 相同问题的解决方案。假设 x 方向有 20 个样本,y 方向有 30 个样本,以及每个交叉点的数据(总共 20x30 = 600 个样本)。我使用x,y = np.mesgrid(x,y),然后使用np.stack((x,y), axis=2).reshape(-1, 2) 得到一个(600,2) 数组,它将是xdata,包含所有600 个x 和y 组合。然后我将 600 个样本中的数据展平为 1d (600,) 数组,该数组将是 ydata 而不是 2d (20, 30) 数组。然后您可以使用x, y = np.hsplit(X, 2) 将您的数据解压缩到func 中。 (Xxdata
    【解决方案2】:

    拟合未知数量的参数

    在这个例子中,我们尝试重现一些测量数据measData。 在此示例中,measData 由函数 measuredData(x, a=.2, b=-2, c=-.8, d=.1) 生成。我练习,我们可能以某种方式测量了measData - 所以我们不知道它是如何在数学上描述的。因此适合。

    我们通过函数polynomFit(inp, *args) 描述的多项式进行拟合。由于我们想尝试不同阶的多项式,因此在输入参数的数量上保持灵活很重要。 自变量(在您的情况下为 x 和 y)在inp 的“列”/第二维中编码。

    import numpy as np
    import matplotlib
    import matplotlib.pyplot as plt
    from scipy.optimize import curve_fit
    
    def measuredData(inp, a=.2, b=-2, c=-.8, d=.1):
        x=inp[:,0]
        y=inp[:,1]
        return a+b*x+c*x**2+d*x**3 +y
    
    def polynomFit(inp, *args):
        x=inp[:,0]
        y=inp[:,1]
        res=0
        for order in range(len(args)):
            print(14,order,args[order],x)
            res+=args[order] * x**order
        return res +y
    
    
    inpData=np.linspace(0,10,20).reshape(-1,2)
    inpDataStr=['({:.1f},{:.1f})'.format(a,b) for a,b in inpData]
    measData=measuredData(inpData)
    fig, ax = plt.subplots()
    ax.plot(np.arange(inpData.shape[0]), measData, label='measuered', marker='o', linestyle='none' )
    
    for order in range(5):
        print(27,inpData)
        print(28,measData)
        popt, pcov = curve_fit(polynomFit, xdata=inpData, ydata=measData, p0=[0]*(order+1) )
        fitData=polynomFit(inpData,*popt)
        ax.plot(np.arange(inpData.shape[0]), fitData, label='polyn. fit, order '+str(order), linestyle='--' )
        ax.legend( loc='upper left', bbox_to_anchor=(1.05, 1))
        print(order, popt)
    
    ax.set_xticklabels(inpDataStr, rotation=90)
    

    结果:

    【讨论】:

      【解决方案3】:

      是的。我们可以为 curve_fit 传递多个变量。我写了一段代码:

      import numpy as np
      x = np.random.randn(2,100)
      w = np.array([1.5,0.5]).reshape(1,2)
      esp = np.random.randn(1,100)
      y = np.dot(w,x)+esp
      y = y.reshape(100,)
      

      在上面的代码中,我生成了 x 一个形状为 (2,100) 的二维数据集,即有两个具有 100 个数据点的变量。我已将因变量 y 与自变量 x 拟合,并带有一些噪音。

      def model_func(x,w1,w2,b):
        w = np.array([w1,w2]).reshape(1,2)
        b = np.array([b]).reshape(1,1)
        y_p = np.dot(w,x)+b
        return y_p.reshape(100,)
      

      我们定义了一个模型函数来建立 yx 之间的关系。
      注意:模型输出的形状函数或预测的y应该是(x的长度,)

      popt, pcov = curve_fit(model_func,x,y)
      

      popt 是一个包含预测参数的一维 numpy 数组。在我们的例子中,有 3 个参数。

      【讨论】:

        【解决方案4】:

        优化具有多个输入维度和可变参数数量的函数

        这个例子展示了如何通过越来越多的系数来拟合具有二维输入 (R^2 -> R) 的多项式。该设计非常灵活,可以为任意数量的非关键字参数定义一次来自curve_fit 的可调用 f。

        最小的可重现示例

        import numpy as np
        import matplotlib.pyplot as plt
        from scipy.optimize import curve_fit
        
        def poly2d(xy, *coefficients):
            x = xy[:, 0]
            y = xy[:, 1]
            proj = x + y
            res = 0
            for order, coef in enumerate(coefficients):
                res += coef * proj ** order
            return res
        
        nx = 31
        ny = 21
        
        range_x = [-1.5, 1.5]
        range_y = [-1, 1]
        target_coefficients = (3, 0, -19, 7)
        
        xs = np.linspace(*range_x, nx)
        ys = np.linspace(*range_y, ny)
        im_x, im_y = np.meshgrid(xs, ys)
        xdata = np.c_[im_x.flatten(), im_y.flatten()]
        im_target = poly2d(xdata, *target_coefficients).reshape(ny, nx)
        
        fig, axs = plt.subplots(2, 3, figsize=(29.7, 21))
        axs = axs.flatten()
        
        ax = axs[0]
        ax.set_title('Unknown polynomial P(x+y)\n[secret coefficients: ' + str(target_coefficients) + ']')
        sm = ax.imshow(
            im_target,
            cmap = plt.get_cmap('coolwarm'),
            origin='lower'
            )
        fig.colorbar(sm, ax=ax)
        
        for order in range(5):
            ydata=im_target.flatten()
            popt, pcov = curve_fit(poly2d, xdata=xdata, ydata=ydata, p0=[0]*(order+1) )
        
            im_fit = poly2d(xdata, *popt).reshape(ny, nx)
        
            ax = axs[1+order]
            title = 'Fit O({:d}):'.format(order)
            for o, p in enumerate(popt):
                if o%2 == 0:
                    title += '\n'
                if o == 0:
                    title += ' {:=-{w}.1f} (x+y)^{:d}'.format(p, o, w=int(np.log10(max(abs(p), 1))) + 5)
                else:
                    title += ' {:=+{w}.1f} (x+y)^{:d}'.format(p, o, w=int(np.log10(max(abs(p), 1))) + 5)
            title += '\nrms: {:.1f}'.format( np.mean((im_fit-im_target)**2)**.5 )
            ax.set_title(title)
            sm = ax.imshow(
                im_fit,
                cmap = plt.get_cmap('coolwarm'),
                origin='lower'
                )
            fig.colorbar(sm, ax=ax)
        
        for ax in axs.flatten():
            ax.set_xlabel('x')
            ax.set_ylabel('y')
        
        plt.show()
        

        附:此答案的概念与我在此处的其他答案相同,但代码示例更加清晰。在给定的时间,我将删除另一个答案。

        【讨论】:

          【解决方案5】:

          是的,有:只需给curve_fit 一个多维数组xData

          【讨论】:

          • 我一直在尝试将 x 和 y 合并为一个数组 z=[x,y], so that x=z[0] 和 y=z[1]。但是curve_fit 似乎不喜欢这样并给了我一个错误:TypeError: unsupported operand type(s) for /: 'list' and 'float'
          猜你喜欢
          • 2020-06-17
          • 2018-03-28
          • 2016-01-27
          • 2017-01-31
          • 2023-03-05
          • 1970-01-01
          • 1970-01-01
          • 1970-01-01
          • 1970-01-01
          相关资源
          最近更新 更多