多项式拟合通过 1 个点，力导数 = 0

Question

我正在寻找类似于numpy.polyfit 的东西 我有一条带有固定点的曲线，看起来像 2 次多项式

我想做的是：

• 准确地通过第一点（在以下示例中为(0.05 , 1.0)）
• 关于第一点有derivative =0

例子：

TabX：

0,050 ; 0,055 ; 0,060 ; 0,065 ; 0,070 ; 0,075 ; 0,080 ; 0,085 ; 0,090 ; 0,095 ;     0,100 ; 0,110 ; 0,120 ; 0,130 ; 0,140 ; 0,150 ; 0,160 ; 0,170 ; 0,180 ; 0,190 ; 0,200 ; 0,210 ; 0,220 ; 0,230 ; 0,243

TabY：

1,000000000 ; 1,000446069 ; 1,000395689 ; 1,000359466 ; 1,000313867 ; 0,999937484 ; 0,999760969 ; 0,999811533 ; 0,999966352 ; 0,999767956 ; 1,000148567 ; 1,000634904 ; 1,000735849 ; 1,001199937 ; 1,001510678 ; 1,001722496 ; 1,001992602 ; 1,002487029 ; 1,003492247 ; 1,004006533 ; 1,004832845 ; 1,005730132 ; 1,006327527 ; 1,007109894 ; 1,008266254

我已经找到了一个“简单但野蛮”的解决方案来解决第一点：我在这一点上添加了很多权重，或者使用权重函数，或者通过在 TabX 中添加大量 0.05 和TabY 中有很多 1.0 并使用正常的 np.polyfit 函数。它很丑，但很管用。

但我真的不知道如何在(0.05;1.0) 点中有derivative=0

另外，在this thread 中，据说拉格朗日乘数可以解决问题，但我无法使用@Jaime 编写的函数，我遇到了错误LinAlgError, 'Singular matrix'

另外，由于我必须能够在基本Abaqus 中启动此脚本，因此该解决方案只能使用基本python 2.7 和numpy。在此实现中无法使用scipy 或matplotlib。

_____________________________________________________________________

编辑：使用 DanielF 的答案，我可以做一些事半功倍的事情（也感谢 DanielF 对我的原始帖子的更正，更易于阅读）

使用新的来源是个好主意 但是，我无法使用我的数据使其高效 这是有效的：

def WorkingDeg4 ():
    x = np.arange(100)
    y0 = 4.0*x**4+0.07 * x ** 3 + 0.3 * x ** 2 + 1.1 * x
    y = y0 + 1000 * np.random.randn(x.shape[0])
    XX = np.vstack((x**4,x ** 3, x ** 2, x, np.ones_like(x))).T
    p_all = np.linalg.lstsq(XX, y)[0]
    pp = np.polyfit(x, y, 3)
    p_no_offset = np.linalg.lstsq(XX[:, :-1], y)[0]
    y_fit = np.dot(p_no_offset, XX[:, :-1].T)
    for i in range(0,len(x)):
        print x[i],y0[i],y[i],y_fit[i]

但是，如果我想使用我的数据制作 int 我主要说的是：

if __name__ == '__main__':
    MyDataX=[0.050 , 0.055 , 0.060 , [...], 0.243]
    MyDataY=[1.000000000 , 1.000446069 , 1.000395689 , [...] , 1.008266254]
    TabX=[0.0]*len(MyDataX)
    TabY=[0.0]*len(MyDataY)
    for i in range(0,len(TabX)):
        TabX[i]=MyDataX[i]-MyDataX[0]
        TabY[i]=MyDataY[i]-MyDataY[0]

所以，在这一点上，我做了“回到原点”阶段

我想和def Working 做同样的事情，但对于我的数据，所以我复制了 tge WorkingDeg4 的过去，只是摆脱了 x 和 y 的创建并将其放入参数中

def NOTWorkingDeg4 (x,y): 
    XX = np.vstack((x**4,x ** 3, x ** 2, x, np.ones_like(x))).T
    p_all = np.linalg.lstsq(XX, y)[0]
    pp = np.polyfit(x, y, 3)
    p_no_offset = np.linalg.lstsq(XX[:, :-1], y)[0]
    y_fit = np.dot(p_no_offset, XX[:, :-1].T)

这个不工作....我有一个神秘的线

XX = np.vstack((x**4,x ** 3, x ** 2, x, np.ones_like(x))).T

说TypeError: unsupported operand type(s) for ** or pow(): 'list' and 'int'

所以据我了解，它不想这样做，因为当它执行 x**4 时，x 不是整数。但是我不知道热解决的问题

_____________________________________________________________________ 编辑2：找到它： 问题是启动 TabX 和 TabY，不是作为数组，而是作为 np.array 错误：

TabX=[0.0]*len(MyDataX)
TabY=[0.0]*len(MyDataY)

正确：

TabX=np.array([0.0]*len(LongueursFissureGlobale))
TabY=np.array([0.0]*len(CourbeInterpolationGlobale))

Answer 1

您可能想要转换您的数据，以便 (0.05, 1.0) 是您的原点 (0, 0) 以使您的计算更容易，然后执行 numpy.linalg.lstsq

x = TabX - 0.05
y = TabY - 1.0
X_poly = np.vstack((x ** 4, x ** 3, x ** 2))
poly_coeffs = np.linalg.lstsq(X_poly.T, y)
y_fit = np.dot(poly_coeffs, X_poly)

如果需要，您必须将多项式转换回您的旧坐标，但这种转换使拟合变得更加简单。

更多详情请见This answer。