最小二乘法理论最小二乘法是一种数学优化技术。利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小。是解决曲线拟合最常用的方法,其思路如下:其中,是预选定的一组线性相关的函数,是待定系数,拟合准则是使与的距离的平方和最小,称为最小二乘法准则。三维空间中,定义4x4的网格,首先定义z值,而网格点上的z值不一样,我们所要做的就是根据这个z值去拟#合这个面上所有点的值。

- 最小二乘法理论
最小二乘法(又称最小平方法)是一种数学优化技术。它通过最小化误差的平方和寻找数据的最佳函数匹配。利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小。是解决曲线拟合最常用的方法,其思路如下:
其中,是预选定的一组线性相关的函数,是待定系数,拟合准则是使与的距离的平方和最小,称为最小二乘法准则。
最小二乘准则进行最小二乘平差计算的一个基本原则
- 代码示例
import pandas as pdimport numpy as npfrom sklearn.linear_model import LinearRegressionimport matplotlib.pyplot as pltfrom mpl_toolkits.mplot3d import Axes3D#----------------------------------------------------------------------------------------------------------------------# Step1:创建需要被拟合的目标。三维空间中,定义4x4的网格,首先定义z值,而网格点上的z值不一样,我们所要做的就是根据这个z值去拟#合这个面上所有点的值。#----------------------------------------------------------------------------------------------------------------------np.random.seed(0)dim = 4Z = (np.ones((dim, dim)) * np.arange(1, dim 1, 1))**3np.random.rand(dim, dim) * 200x = np.arange(1, dim 1).reshape(-1, 1)y = np.arange(1, dim 1).reshape(1, -1)X, Y = np.meshgrid(x, y)#----------------------------------------------------------------------------------------------------------------------# Step2:自定义一组线性相关的函数, 3阶。#----------------------------------------------------------------------------------------------------------------------features = {}features['x^0*y^0'] = np.matmul(x**0, y**0).flatten()features['x*y'] = np.matmul(x, y).flatten()features['x*y^2'] = np.matmul(x, y**2).flatten()features['x^2*y^0'] = np.matmul(x**2, y**0).flatten()features['x^2*y'] = np.matmul(x**2, y).flatten()features['x^3*y^2'] = np.matmul(x**3, y**2).flatten()features['x^3*y'] = np.matmul(x**3, y).flatten()features['x^0*y^3'] = np.matmul(x**0, y**3).flatten()dataset = pd.DataFrame(features)#----------------------------------------------------------------------------------------------------------------------# Step3:将选定函数与目标值带入SkLearn包中的线性回归拟合模块,它可以使平方和最小,结果返回截距和斜率。#----------------------------------------------------------------------------------------------------------------------reg = LinearRegression().fit(dataset.values, Z.flatten())# reg.intercept_为截距, reg.coef_为斜率z_pred = reg.intercept_np.matmul(dataset.values, reg.coef_.reshape(-1, 1)).reshape(dim, dim)#----------------------------------------------------------------------------------------------------------------------# Step4:可视化。#----------------------------------------------------------------------------------------------------------------------fig = plt.figure(figsize=(5, 5))ax = Axes3D(fig)ax.plot_surface(X, Y, z_pred, label='prediction', cmap=plt.get_cmap('rainbow'))ax.scatter(X, Y, Z, c='r', label='datapoints')plt.show()
结果如下:
图1
上例定义的多项式阶数为3,对于大多数问题已经足够了,如果想定义更高阶数,则可参考如下代码:
import itertoolsimport numpy as npimport matplotlib.pyplot as pltfrom mpl_toolkits.mplot3d import Axes3D#-----------------------------------------------------------------------------------------------------------------------# Step1: 创建线性相关的函数,阶数自己定义,结果为拟合系数;#-----------------------------------------------------------------------------------------------------------------------def polyfit2d(x, y, z, order):ncols = (order1)**2G = np.zeros((x.size, ncols))ij = itertools.product(range(order 1), range(order 1))for k, (i, j) in enumerate(ij):G[:, k] = x**i * y**jm, _, _, _ = np.linalg.lstsq(G, z, rcond=-1) # lstsq的输出包括四部分:回归系数、残差平方和、自变量X的秩、X的奇异值return m#-----------------------------------------------------------------------------------------------------------------------# Step2: 创建拟合函数,将欲拟合值和拟合系数待入,返回预测值;#-----------------------------------------------------------------------------------------------------------------------def polyval2d(x, y, m):order = int(np.sqrt(len(m))) - 1 # 根据多项式的列数反算阶数ij = itertools.product(range(order 1), range(order 1))z = np.zeros_like(x)for a, (i, j) in zip(m, ij):z= a * x**i * y**jreturn z#-----------------------------------------------------------------------------------------------------------------------# Step3: 示例;#-----------------------------------------------------------------------------------------------------------------------x = np.array([4, 5, 5, 4])y = np.array([2, 3, 4, 5])z = np.array([2, 3, 4, 7])N_ORDER = 4m = polyfit2d(x, y, z, N_ORDER)N_MESH = 10xx, yy = np.meshgrid( np.linspace(x.min(), x.max(), N_MESH),np.linspace(y.min(), y.max(), N_MESH))zz = polyval2d(xx, yy, m)#-----------------------------------------------------------------------------------------------------------------------# Step4: 可视化;#-----------------------------------------------------------------------------------------------------------------------fig = plt.figure(figsize=(5, 5))ax = Axes3D(fig)ax.plot_surface(xx, yy, zz, label='prediction', cmap=plt.get_cmap('rainbow'))ax.scatter(x, y, z, c='r', label='datapoints')plt.show()
结果如下:
最小二乘法很基本也很常用,其本质是插值,但是我发现在当数值非常大时它的效果却不是很好,这时候就需要用到一些其他的方法,比如克里金法等。
声明:仅供参考
