我们尝试实现一个标准的多项式拟合程序，期间遇到一个无法理解的问题。下面是我们的代码，其中包含了样本数据 x 和 y，我们使用正规方程和 numpy 中的 polyfit 函数对数据进行拟合。
def mytest(N):
print "\nPolynomial order (N): {}".format(N)
mVals = [5, 10, 15, 20]
a_orig = [198.764, 13.5, 0.523]
for M in mVals:
x = arange(-M, M+1)
y = matrix(a_orig[0]+ a_orig[1]*x + a_orig[2]*x**2).T
# Code implementing the solution from the normal equations
nArray = arange(N+1)
A = matrix([[n**i for i in nArray] for n in x])
B = (A.T*A).I
a_myfit = B*(A.T*y)
# numpy's polyfit
a_polyfit = polyfit(x, y, N)
print "M: {}".format(M)
print ["{0:0.3f}".format(i) for i in a_orig]
print ["{0:0.3f}".format(i) for i in array(a_myfit)[:,0]]
print ["{0:0.3f}".format(i) for i in list(array(a_polyfit)[:,0])[::-1]]
mytest(N=5)
mytest(N=6)
以下是我们尝试使用 N=5 和 N=6 时得到的输出结果：
Polynomial order (N): 5
M: 5
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000']
['198.764', '13.500', '0.523', '-0.000', '0.000', '0.000']
M: 10
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '-0.000', '0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '0.000']
M: 15
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000']
['198.764', '13.500', '0.523', '-0.000', '-0.000', '0.000']
M: 20
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '-0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '0.000']
Polynomial order (N): 6
M: 5
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000', '-0.000']
['198.764', '13.500', '0.523', '-0.000', '0.000', '0.000', '-0.000']
M: 10
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '-0.000', '-0.000', '-0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000', '-0.000']
M: 15
['198.764', '13.500', '0.523']
['294.451', '13.500', '-0.061', '0.000', '-0.001', '-0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '0.000', '-0.000']
M: 20
['198.764', '13.500', '0.523']
['369.135', '13.500', '-0.046', '0.000', '-0.000', '-0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000', '-0.000']
从输出中可以看到，当 N > 5 且 M > 13 时，多项式拟合的结果出现了错误。
2、解决方案
为了解决这个问题，我们需要检查代码中的问题并进行修改。主要需要注意以下几点：
- 确保矩阵 A 的数据类型能够处理较大的数值，以避免出现溢出。
- 检查数据类型是否一致，特别是当使用多个库或模块时，确保数据类型之间兼容。
通过检查代码，发现问题出现在矩阵 A 的数据类型上。在默认情况下，A 的数据类型是 int32，它只能表示有限范围的整数。当指数非常大时，n**i 的值可能会超出 int32 的表示范围，导致溢出并产生不正确的结果。
要解决这个问题，我们可以将 A 的数据类型更改为 float64，它能够表示更大的数值。修改后的代码如下：
A = np.matrix([[n**i for i in nArray] for n in x], dtype='float64')
修改后，我们使用相同的数据重新运行代码，得到了正确的结果：
Polynomial order (N): 5
M: 5
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '0.000', '-0.000']
['198.764', '13.500', '0.523', '-0.000', '0.000', '0.000']
M: 10
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '-0.000', '0.000', '-0.000']
['198.764', '13.500', '0.523', '0.000', '0.000', '0.000']
M: 15
['198.764', '13.500', '0.523']
['198.764', '13.500', '0.523', '0.000', '0.0
