[数值计算-17]：最小二乘法的求解2- 链式求导与多元线性方程组求解

Posted 2021-09-12 文火冰糖的硅基工坊

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了[数值计算-17]：最小二乘法的求解2- 链式求导与多元线性方程组求解相关的知识，希望对你有一定的参考价值。

作者主页(文火冰糖的硅基工坊)：https://blog.csdn.net/HiWangWenBing

本文网址：https://blog.csdn.net/HiWangWenBing/article/details/119977443

3.2 线性拟合函数的最小二乘法的几何意义

3.3 最小二乘算法实现

3.4 最小二乘法求解1：“法向量”方程组法

3.5 最小二乘法求解1：二元线性方程组示例

第4章最小二乘法求解2：多元线性方程组解法

4.0 前置条件：

4.1 步骤1：构建样本数据集

4.2 步骤2：构建拟合函数

4.3 步骤3：构建最小二乘的残差函数，并计算最佳参数（最关键步骤）

4.4 步骤4：利用拟合函数进行数据预测

4.5 步骤5：图形展示

第1章什么是函数逼近？

1.1 抛出问题1：函数插值

利用有限的样本数据，发现其内在的规律，并用这个规律预测未来新的数据。

（1）单个数据点

0次函数通过样本点：唯一确定一个点 y = f(x) = a0
1次直线函数通过样本点：可以有无数 y = f(x) = a1x + a0
2次抛物线函数通过样本点：可以有无数 y = f(x) = a1x^2 + a1x + a0

（2）2个数据点

0次函数通过样本点：无
1次直线函数通过样本点：唯一直线 y = f(x) = a1x + a0
2次抛物线函数通过样本点：可以有无数 y = f(x) = a1x^2 + a1x + a0

（3）3个数据点

0次函数通过样本点：无
1次直线函数通过样本点：无
2次抛物线函数通过样本点：唯一抛物线y = f(x) = a1x^2 + a1x + a0

问题：

如果有(xn+1,yn+1), (xn,yn)........(x1,y1), (x0,y0)样本点，那么如何选择一个最低次的多项式函数，可以穿过上述样本点？

推测：

对应n+1个点，可以唯一的确定一个一元n次的多项式函数，该多项式函数可以穿越所有n+1个点。

1.2 抛出问题2：函数拟合

如果有n个点，不要求选出的函数穿越所有的点，而是根据这些点构建的轮廓，选择一个更低维度（次数）的函数尽可能的靠近这些样本点呢？

当函数的次数远远小于样本点的次数是，该如何选择低次的函数?

1.3 什么是回归

回归与拟合是基本相同的概念。线性拟合通常称为线性回归。

第2章函数拟合

2.1 插值有时候插值并不是最好的选择

（1）简单直线拟合好于多项式拟合

在上图图有，有几十个样本点，如果用多项式插值来进行函数逼近，虽然插值函数在样本点的误差为0，但很显然，需要高次多项式来进行逼近。

如果样本点有几百个，那么多项式的次数高达几百次。

从图中可以大致看出，要逼近上述样本点，只需要一个简单的线性函数即可：y=ax+b，这种函数逼近称为拟合。

（2）抛物线逼近

同理，对上述样本点的逼近，也不需要几十次的多项式进行插值，值需要一个简单的抛物线就可以进行逼近：y =ax^2 + bx + c，，这种函数逼近称为拟合。

（3）三角函数拟合

2.2 什么函数拟合

2.2 函数拟合的基本步骤

步骤1：构建样本

（1）读取已有样本数据，并图形展示样本数据。或者

（2）构建自行测试样本数据

（3）图形化展示样布数据（散点图）

步骤2：选择函数模型

根据散点图，选择函数类，函数类可以从初等函数中进行选取，如线性函数、二次或多次多项式函数、三角函数等。

步骤3：构建最小二乘的残差函数，并计算最佳参数（最关键步骤）

（1）选择最佳拟合的范数（这里选择最小二乘法）

选取了函数类型后，每个函数都有自己的待定参数，不同的参数，其拟合效果是不同的，如何选择函数的参数，使得拟合效果最好呢？常有的有最小二乘法。

最小二乘法是函数拟合最重要的环节之一。

（2）计算最佳参数（这是最关键、最核心的步骤）

最小二乘法得到的loss残差函数是一个多元二次多项式函数：

残差函数的“元”的个数，就是拟合函数的参数
残差函数的“元”的次数：2次
拟合最好的依据：多元函数的最小值点处的各个变量值，就是拟合函数最好的参数。

求解多元二次多项式函数的方法有多种：

法向量解线性方程组法
梯度下降法

无论哪种方法，都是通过求残差函数的最小值，基于指定样本数据，获得了最佳参数的拟合函数。

步骤4：利用拟合函数进行数据预测

生成一组输入数据，利用第4步获得的拟合函数的参数以及对应的拟合函数，对数据进行函数运算（预测），得到新的输出数据。

步骤5：图形展示

把预测数据与采样数据、理论曲线一起进行图形展示。

第3章最小二乘法的线性方程组解法

3.1 什么是最小二乘法

3.2 线性拟合函数的最小二乘法的几何意义

注意：

最小二乘选择一个适当的直线参数，使得所有样本点处的yi - f(xi) = yi -（axi+b)的平方和的平均值最小，而不是所有样本点到拟合直线的距离的和的平均值最小。

这样选择，主要处于如下原因：

（1）计算简单：求点到拟合曲线的距离要比直接计算样本点处的Y值的差复杂很多。

（2）符合拟合函数原本的物理意义：拟合函数要在所有的样本点xi处，其f(xi）也尽可能的相等，或者说误差尽可能的小。

3.3 最小二乘算法实现

（1）最小二乘法求解1：一元参数的解析法求解

（2）最小二乘法求解2：“法向量”线性方程组求解法

（3）最小二乘法求解3：梯度下降法

最小二乘（二范数）的优点：

一阶和二阶导数的存在，使得该残差函数，很容易通过“梯度下降法”获得其最小值和对应的各个参数的值。

3.4 最小二乘法求解1：“法向量”方程组法

3.5 最小二乘法求解1：二元线性方程组示例

第4章最小二乘法求解2：多元线性方程组解法

4.0 前置条件：

#导入库
from math import *
import time
import numpy as np
import matplotlib.pyplot as plt #画图工具
from pylab import mpl           #中文字体
from scipy import optimize      #最小二乘算法的算法库

4.1 步骤1：构建样本数据集

#步骤1：构建样本

#(1) 采用np, 直接手工生成样本的输入：一组等距离的分布在[-1,1]之间的100个点
sample_numbers = 50

x_data = np.linspace(-1, 1, sample_numbers)

#(2) 为这些数据手工打上理论输出值（标签值）：y = 2x + 1
y_data_pure =  2 * x_data + 1.0

#（3）为了模拟现实情况，通过随机数来模拟数据噪声
noise_range = 0.4
np.random.seed(10) #设置随机种子, 确保不同时候，执行结果是相同的
#randn(n)生成的0为均值，1为标准差的正态分布的n个随机数。
y_noise = np.random.randn(*x_data.shape) * noise_range  # *x_data.shape：输入样本的维度或个数

#（4）人工生成样本的输出：理论值 + 噪声
y_data_noise = y_data_pure + y_noise

#(5) 显示样本数据
# 样本的散点图
plt.scatter(x_data, y_data_noise, label="sample", color="black")

# 内在的、理论的曲线图
plt.plot(x_data, y_data_pure, label="f_pure(x)", color="blue", linewidth = 4)

#设置属性
mpl.rcParams['font.sans-serif'] = ['SimHei']
mpl.rcParams['axes.unicode_minus'] = False
plt.title("线性拟合")
plt.legend(loc="upper left")
plt.show()

4.2 步骤2：构建拟合函数

#步骤2：构建拟合函数：二元一次拟合函数
def f(x, a, b):
    return (a*x + b)

4.3 步骤3：构建最小二乘的残差函数，并计算最佳参数（最关键步骤）

#步骤3：利用python库提供的最小二乘算法来计算拟合函数的参数
rst_fit = optimize.curve_fit(f, x_data, y_data_noise)
print(rst_fit)

a_fit, b_fit = rst_fit[0]
print("参数a=", a_fit)
print("参数b=", b_fit)

(array([1.95913373, 1.04099449]), array([[7.76417539e-03, 1.62852841e-11],
       [1.62852841e-11, 2.69369357e-03]]))
参数a= 1.9591337280301642
参数b= 1.0409944876382622

备注说明：

这里利用python库提供的最小二乘算法来计算拟合函数的参数。

因此没有展现用数值的方法求最小二乘最佳参数的过程。

如下函数就展现通过求解多元线性方程组的方法，获得拟合函数的参数a和b。

#自定义最小二乘求解拟合函数参数：偏导+线性方程组（求偏导后一次函数）
def usr_curve_fit(f, x_data, y_data_noise):
    #获取采样数据的个数
    n = len(x_data)
    
    #计算线性方程组矩阵的数据
    sum_xi = np.sum(x_data)
    sum_yi = np.sum(y_data_noise)
    sum_xiyi = np.sum(x_data * y_data_noise)
    sum_xixi = np.sum(x_data**2)
    
    #构建线性方程组矩阵
    A=[[0,0],[0,0]]
    b=[0,0]
    A[0][0] = n
    A[0][1] = sum_xi
    A[1][0] = sum_xi
    A[1][1] = sum_xixi
    b[0] = sum_yi
    b[1] = sum_xiyi
    
    #解方程组
    rst = np.linalg.solve(A,b)
    
    #返回拟合函数参数
    return (rst[1],rst[0])

rst_usr = usr_curve_fit(f, x_data, y_data_noise)
print(rst_usr)

a_usr, b_usr = rst_usr
print("参数a=", a_usr)
print("参数b=", b_usr)

输出：
(1.959133728205619, 1.040994488066939)
参数a= 1.959133728205619
参数b= 1.040994488066939

可以看出，通过自定义的线性方程组求解获取的拟合参数的参数与通过scipy获得的参数值是一致的。

参数a= 1.9591337280301642 参数b= 1.0409944876382622

参数a= 1.959133728205619 参数b= 1.040994488066939

4.4 步骤4：利用拟合函数进行数据预测

# 步骤4：利用获得的拟合函数进行数据预测

# scipy算法的拟合数据
y_data_fit = f(x_data, a_fit, b_fit)

# 线性方程组求解的拟合数据
y_data_usr = f(x_data, a_usr, b_usr)

4.5 步骤5：图形展示

#步骤5： 图形化展示
#(1) 显示样本数据曲线
plt.scatter(x_data, y_data_noise, label="sample", color="black")

#(2) 显示理论数据曲线
plt.plot(x_data, y_data_pure, label="intrinsic", color="blue", linewidth = 2)

#(3-1) 显示预测数据曲线 -  scipy库实现
plt.plot(x_data, y_data_fit,  label="predict",  color="red",  linewidth = 2)

#(3-2) 显示预测数据曲线 -  自定义实现
plt.plot(x_data, y_data_usr,  label="predict",  color="green",  linewidth = 2)

#设置属性
mpl.rcParams['font.sans-serif'] = ['SimHei']
mpl.rcParams['axes.unicode_minus'] = False
plt.title("线性拟合")
plt.legend(loc="upper left")
plt.show()