带适配参数的numpy.polyfit

Question

关于这个问题：多项式方程的参数，我得到了一个平方函数的三个参数，公式是 y = a*x² + b*x + c。现在我只想得到这个平方函数的第一个参数，也就是描述我的函数 y = a*x² 的参数。换句话说，我想把 b=c=0，然后得到适合的 a 参数。如果我理解得没错，polyfit 这个工具是无法做到这一点的。

Answer 1

这些系数是用来最小化平方误差的，你并不是自己去指定它们的值。不过，如果某些系数的影响太小，你可以把它们设为零。比如，我有一组在曲线 y = 33*x² 上的点：

In [51]: x=np.arange(20)

In [52]: y=33*x**2  #y = 33*x²

In [53]: coeffs=np.polyfit(x, y, 2)

In [54]: coeffs
Out[54]: array([  3.30000000e+01,   8.99625199e-14,  -7.62430619e-13])

In [55]: epsilon=np.finfo(np.float32).eps

In [56]: coeffs[np.abs(coeffs)<epsilon]=0

In [57]: coeffs
Out[57]: array([ 33.,   0.,   0.])

Answer 2

这可以通过 numpy.linalg.lstsq 来实现。为了说明如何使用它，可能最简单的方式是先手动演示一下如何进行标准的二次多项式拟合。假设你有测量数据向量 x 和 y，首先你需要构建一个叫做 设计矩阵 M，其构造方式如下：

M = np.column_stack((x**2, x, np.ones_like(x)))

接下来，你可以通过最小二乘法来求解方程 M * k = y，使用 lstsq 方法，具体如下：

k, _, _, _ = np.linalg.lstsq(M, y)

这里 k 是一个列向量 [a, b, c]，包含了通常的系数。需要注意的是，lstsq 会返回一些其他参数，你可以选择忽略它们。这是一个非常强大的技巧，它可以让你将 y 拟合到你放入设计矩阵的列的任意线性组合中。比如，它可以用于类型为 z = a * x + b * y 的二维拟合（例如，这个例子，我在 Matlab 中使用了同样的技巧），或者处理像你问题中那样缺失系数的多项式拟合。

在你的例子中，设计矩阵其实就是一个包含 x**2 的单列。快速示例：

import numpy as np
import matplotlib.pylab as plt

# generate some noisy data
x = np.arange(1000)
y = 0.0001234 * x**2 + 3*np.random.randn(len(x))

# do fit
M = np.column_stack((x**2,)) # construct design matrix
k, _, _, _ = np.linalg.lstsq(M, y) # least-square fit of M * k = y

# quick plot
plt.plot(x, y, '.', x, k*x**2, 'r', linewidth=3)
plt.legend(('measurement', 'fit'), loc=2)
plt.title('best fit: y = {:.8f} * x**2'.format(k[0]))
plt.show()

结果：