拟合函数未知数个数与用于拟合的序列点数的关系

• 1.以二元一次函数为例（a、b为未知数），当序列点数为1个点或2个点，即序列点数小于拟合函数未知数个数或大于等于拟合函数未知数个数时，用最小二乘法构造目标函数，求其驻点，得到最优的拟合函数系数：

从以上推导可知：

• 当拟合函数有两个系数，而只有一个点可用于拟合时，虽然在求解最优系数时，仍有两个偏导数方程，但这两个方程线性相关，其实只有一个方程是有效的，因为只有一个点进行拟合，相当于一个方程无法求解两个未知数，所以这时系数a、b有无穷多解；表现在几何上则是过一点有无穷多条最优直线，使真实值与拟合值的差值平方和最小。

• 当有两个点用于拟合时，这两个方程线性无关，系数矩阵行满秩，相当于两个方程可以求解两个未知数，这时a、b有唯一解；表现在几何上则是过两点有唯一一条最优直线，使真实值与拟合值的差值平方和最小。

• 当有多个点用于拟合时，这两个方程线性无关，系数矩阵行满秩，相当于两个方程可以求解两个未知数，这时a、b有唯一解；表现在几何上则是过多点有唯一一条最优直线，使真实值与拟合值的差值平方和最小。

• 2.下面以python中多项式拟合为例，使拟合函数系数的个数分别小于或大于用于拟合的序列点数，查看各自的拟合图形和返回参数，并作出说明。

import matplotlib.pyplot as plt
import numpy as np
import copy
plt.style.use('seaborn-white')


x = [np.linspace(0, 35, 36), np.linspace(0, 9, 10), np.linspace(0, 19, 20)]
# 使绘图时的x取值比多项式拟合时，自变量x的取值更密集，以观察到可能出现的振荡
x_plot = [np.linspace(0, 35, 36*5), np.linspace(0, 9, 10*5), np.linspace(0, 19, 20*5)]
y = [np.random.rand(len(x[0])), np.random.randn(len(x[1]))*10, np.random.randn(len(x[2]))*10]
y[0][int(len(x[0])/2-2):int(len(x[0])/2+2)] = 10
n = [[2, 3, 4], [4, 9, 19], [9, 19, 39]]
weights = [(max(y[i])*2-y[i])/sum(max(y[i])*2-y[i]) for i in range(len(x))]

c = copy.deepcopy(n)  # c与n具有相同的结构；一定要用copy.deepcopy，不能用浅拷贝copy，否则后面对c的赋值也会改变n
p = copy.deepcopy(n)  # p与n具有相同的结构
for i in range(len(x)):
    for j in range(len(n[i])):
        c[i][j] = np.polyfit(x[i], y[i], n[i][j], full=True, w=weights[i])
        print('序列长度{0}个点，拟合次数为{1}，即目标函数未知数个数为{2}，即目标函数对自变量偏导数的方程个数为{2}：'.format(len(x[i]), n[i][j], n[i][j]+1))
        print('拟合函数的系数（次数由高到低）：', '\n', c[i][j][0])
        print('拟合值与实际值的MAPE：', sum(abs((p[i][j](x[i])-y[i])/y[i]))/len(x[i]))
        print('目标函数对自变量偏导数方程组的系数矩阵的秩：', c[i][j][2])
        print('目标函数对自变量偏导数方程组的系数矩阵的奇异值：', '\n', c[i][j][3])
        print('拟合的相关条件数：', c[i][j][4], '\n')
        p[i][j] = np.poly1d(c[i][j][0])
    plt.figure(figsize=(10, 6))
    plt.scatter(x[i], y[i], facecolors='none', edgecolor='darkblue', label='original_data')
    plt.plot(x[i], y[i], color='darkblue')
    plt.plot(x_plot[i], p[i][0](x_plot[i]), '-', color='g', label='polynomial of degree {}'.format(n[i][0]))
    plt.plot(x_plot[i], p[i][1](x_plot[i]), '--', color='orange', label='polynomial of degree {}'.format(n[i][1]))
    plt.plot(x_plot[i], p[i][2](x_plot[i]), color='magenta', label='polynomial of degree {}'.format(n[i][2]))
    plt.ylim(min(y[i])-abs(min(y[i])/2), max(y[i])*1.5)  # 使y坐标的上下限关于y的最大最小值对称，且不受多项式拟合函数的振荡值影响
    plt.legend()
    plt.title('series length: {0}, degree of polyfit: {1}, {2}, {3}'.format(len(x[i]), n[i][0], n[i][1], n[i][2]))


#########################################################################
def polynomial_all(data, n):  # data为df['amou_crct']
    """
    param：
        data：理论销售金额，根据节日效应期长度，向左右各取2倍，共5倍长度，若遇节日则剔除补足，并在data全长内剔除周日效应
        n: 拟合次数
    return：
        fitts：经多项式拟合后的各点拟合值
    """
    if len(data) <= 3:
        raise Exception('历史数据过少，不能计算有效的节日效应系数')
    else:
        Y = np.array(data)
        X = np.linspace(1, len(Y), len(Y))  # 等差数列
        weights = (max(Y) * 2 - Y) / sum(max(Y) * 2 - Y)

        C = np.polyfit(X, Y, n, full=True, w=weights)  # C:多项式拟合的系数，其对应的次数从高到低排列
        p = np.poly1d(C[0])
        fitts = p(X)

    return fitts


x1 = list(range(1, 26))
y1 = np.random.rand(len(x1))
y1[int(len(x1)/2-2):int(len(x1)/2+1)] = 5
n = 3

fitts = polynomial_all(y1, n)

plt.figure('序列长度{0}个点，拟合次数为{1}'.format(len(x1), n))
plt.plot(x1, y1, '.', x1, fitts, '-')
plt.ylim(min(y1)-abs(min(y1)/10), max(y1)+max(y1)/10)

• 3.序列长度10个点，拟合次数为4，即目标函数未知数个数为5，即目标函数对自变量偏导数的方程个数为5：
  
  

• 当序列点数为10，拟合次数为4时，有5个系数，目标函数对自变量偏导数的方程有5个，目标函数对自变量偏导数方程组的系数矩阵的秩为5，系数矩阵行满秩，有5个奇异值；由于行满秩，此时若调用求解方法通常可得到系数的唯一实数解；np.polyfit内部调用的局部优化的梯度下降算法也能得到较合理的系数，从而得到较好的拟合结果。

• 当序列点数为10，拟合次数为19时，有20个系数，目标函数对自变量偏导数的方程有20个，目标函数对自变量偏导数方程组的系数矩阵的秩为10，系数矩阵行不满秩，有10个奇异值；由于行不满秩，此时若调用求解方法无法解出系数，方程组无解或有无穷多解；np.polyfit内部调用的局部优化的梯度下降算法也不能得到合理的系数，从而使拟合区间内部出现剧烈震荡。

• 4.序列长度20个点，拟合次数为9，即目标函数未知数个数为10，即目标函数对自变量偏导数的方程个数为10：
  

• 当序列点数为20，拟合次数为9时，有10个系数，目标函数对自变量偏导数的方程有10个，目标函数对自变量偏导数方程组的系数矩阵的秩为10，系数矩阵行满秩，有10个奇异值；由于行满秩，此时若调用求解方法可得到系数的唯一解；np.polyfit内部调用的局部优化的梯度下降算法也能得到较合理的系数，从而得到较好的拟合结果。

• 当序列点数为20，拟合次数为39时，有40个系数，目标函数对自变量偏导数的方程有40个，目标函数对自变量偏导数方程组的系数矩阵的秩为20，系数矩阵行不满秩，有20个奇异值；由于行不满秩，此时若调用求解方法无法解出系数，方程组无解或有无穷多解；np.polyfit内部调用的局部优化的梯度下降算法也不能得到合理的系数，从而使拟合区间内部出现剧烈震荡。

• 5.总结

• 当未知数个数<<序列点数时，目标函数对自变量偏导数方程组的系数矩阵行满秩，调用求解方法通常能求出方程组的实数解，相当于用n个方程解n个未知数，在拟合时调用的局部优化的梯度下降算法也能得到较好的系数，从而得到较好的拟合效果。

• 当未知数个数＞序列点数时，目标函数对自变量偏导数方程组的系数矩阵行不满秩，调用求解方法不能求出方程组的解，相当于用n个方程不能解出大于n个的未知数，在拟合时调用的局部优化的梯度下降算法也不能得到合理的系数，在拟合区间内部就会出现剧烈震荡。

• 当未知数个数＝序列点数，且未知数个数较少时，目标函数对自变量偏导数方程组的系数矩阵行满秩，调用求解方法通常能求出方程组的实数解，相当于用少数n个方程求解少数n个未知数，在拟合时调用的局部优化的梯度下降算法通常能得到较合理的系数，拟合效果通常较好。

• 当未知数个数＝序列点数，且未知数个数较多时，目标函数对自变量偏导数方程组的系数矩阵行满秩，由于未知数较多，调用求解方法通常不能求出方程组的实数解，在拟合时调用的局部优化的梯度下降算法通常不能得到较合理的系数，拟合效果通常不好，在拟合区间内部易出现剧烈震荡。