MATLAB编写高等数值大作业最小二乘法

最小二乘法

已知有如下数据：

x	4	4.2	4.5	4.7	5.1	5.5	5.9	6.3	6.8	7.1
y	102.56	113.18	130.11	142.05	167.53	195.14	224.87	256.73	299.50	326.72

（1） 分别构造一次、二次、三次最小二乘多项式，画图并计算平方误差；

（2） 构造形如  和的最小二乘多项式，画图并计算平方误差.

1. 数学原理

最小二乘法是一种数学优化技术。它通过最小化误差的平方和寻找数据的最佳函数匹配。利用最小二乘法可以简便地求得未知的数据，并使得这些求得的数据与实际数据之间误差的平方和为最小。最小二乘法还可用于曲线拟合。其他一些优化问题也可通过最小化能量或最大化熵用最小二乘法来表达。

从整体上考虑近似函数 同所给数据点 (i=0,1,…,m)误差 (i=0,1,…,m)的大小，常用的方法有以下三种：一是误差(i=0,1,…,m)绝对值的最大值 ，即误差向量 的∞—范数；二是误差绝对值的和 ，即误差向量r的1—范数；三是误差平方和 的算术平方根，即误差向量r的2—范数；前两种方法简单、自然，但不便于微分运算 ，后一种方法相当于考虑 2—范数的平方，因此在曲线拟合中常采用误差平方和 来度量误差 (i=0，1，…，m)的整体大小。

数据拟合的具体作法是：对给定数据 (i=0,1,…，m)，在取定的函数类 中,求 ,使误差 (i=0,1,…,m)的平方和最小，即

    

从几何意义上讲，就是寻求与给定点 (i=0,1,…,m)的距离平方和为最小的曲线 。函数p(x) 称为拟合函数或最小二乘解，求拟合函数p(x)的方法称为曲线拟合的最小二乘法。

2.2 程序设计

2.2.1（1）问程序代码

% x、 y为输入量；a[num]为拟合函数的系数列向量；x0、y0为拟合函数曲线；delta2为平方误差.

x=[4,4.2,4.5,4.7,5.1,5.5,5.9,6.3,6.8,7.1];

y=[102.56,113.18,130.11,142.05,167.53,195.14,224.87,256.73,299.50,326.72];

[~,k]=size(x);

for n=1:3

    X0=zeros(n+1,k);

    for k0=1:k           %构造矩阵X0

        for n0=1:n+1

            X0(n0,k0)=x(k0)^(n+1-n0);

        end

    end

    X=X0';

    a=(X'*X)\X'*y';

    for i=1:n+1          %answer矩阵存储每次求得的方程系数，按列存储

       answer(i,n)=a(i);

    end

   x0=3.5:0.01:7.5;

    y0=a(1)*x0.^n    ;%根据求得的系数初始化并构造多项式方程

    for num=2:1:n+1    

        y0=y0+a(num)*x0.^(n+1-num);

    end

    subplot(2,2,n)

    plot(x,y,'*')

    hold on

    plot(x0,y0)

      delta2=0;%平方误差

    for i=1;1;10;

        y1(i)=a(1)*x(i).^n;

         for num=2:1:n+1

            y1=y1+a(num)*x(i).^(n+1-num);

         end

    delta2=delta2+(y(i)-y1)^2;

    end

    delta2

end

suptitle('一次、二次、三次最小二乘多项式方程曲线拟合结果')

2.2.2 （2.1）问程序代码

% x、 yy为输入量；y为y=lnyy函数值转换;a、b为拟合函数的系数；x0、y0为拟合函数曲线；delta2为平方误差.

x=[4,4.2,4.5,4.7,5.1,5.5,5.9,6.3,6.8,7.1];

yy=[102.56,113.18,130.11,142.05,167.53,195.14,224.87,256.73,299.50,326.72];

y=log(yy);

e=2.718281828459;

[~,k]=size(x);

X0=zeros(2,k);

    for k0=1:k           %构造矩阵X0

        for n0=1:2

            X0(n0,k0)=x(k0)^(2-n0);

        end

    end

    X=X0';

    a=(X'*X)\X'*y';

    for i=1:2          %answer矩阵存储每次求得的方程系数，按列存储

       answer(i,1)=a(i);

    end

   x0=3.5:0.01:7.5;

    y0=a(1)*x0;   %根据求得的系数初始化并构造多项式方程

    y0=y0+a(2);

    b=e^a(2);%a(2)=lnb;

    a=a(1);

    y=e.^y;

    y0=e.^y0;

    plot(x,yy,'*')

    hold on

    plot(x0,y0)

    delta2=0;%平方误差

    for i=1;1;10;

        y1(i)=b*exp(a*x(i));

    delta2=delta2+(yy(i)-y1)^2;

    end

    delta2

suptitle('构造y=b*exp(a*x)最小二乘多项式')

2.2.2 （2.2）问程序代码

% xx、 yy为输入量；x=lnxx、y=lnyy为变量、函数值转换;a、b为拟合函数的系数；x0、y0为拟合函数曲线；delta2为平方误差.

xx=[4,4.2,4.5,4.7,5.1,5.5,5.9,6.3,6.8,7.1];

yy=[102.56,113.18,130.11,142.05,167.53,195.14,224.87,256.73,299.50,326.72];

y=log(yy);

x=log(xx);

e=2.718281828459;

[~,k]=size(x);

X0=zeros(2,k);

    for k0=1:k           %构造矩阵X0

        for n0=1:2

            X0(n0,k0)=x(k0)^(2-n0);

        end

    end

    X=X0';

    a=(X'*X)\X'*y';

    for i=1:2          %answer矩阵存储每次求得的方程系数，按列存储

       answer(i,1)=a(i);

    end

   x0=x(1):0.01:x(10);

    y0=a(1)*x0;   %根据求得的系数初始化并构造多项式方程

    y0=y0+a(2);

    b=e^a(2);%a(2)=lnb;

    a=a(1);

    x0=e.^x0;

    y0=e.^y0;

    plot(xx,yy,'*')

    hold on

    plot(x0,y0)

    delta2=0;%平方误差

    for i=1;1;10;

        y1(i)=b*xx(i)^a;

    delta2=delta2+(yy(i)-y1)^2;

    end

    delta2

suptitle('构造y=b*x^a最小二乘多项式')

2.3结果分析和讨论

2.3.1下图分别为一次、二次、三次最小二乘多项式方程曲线拟合结果：

 

一次最小二乘多项式：y= 72.0845x -194.1382.

平方误差为69.8924.

二次最小二乘多项式：y= 6.6182x2 -1.1435x+ 1.2356.

平方误差为5.1243e-05.

三次最小二乘多项式：y= -0.0137 x3+ 6.8456 x2 -2.3792x+ 3.4291.

平方误差为3.9466e-05.

 2.3.2下图为y=b*exp(a*x)最小二乘多项式拟合曲线:

 

y=b*exp(a*x)最小二乘多项式：  .

平方误差为25.2762.

2.3.3下图为y=b*x^a最小二乘多项式拟合曲线:

 

y=b*exp(a*x)最小二乘多项式：  .

平方误差为3.2377e-05.

如需要代码请在此链接下下载：https://download.csdn.net/download/weixin_41788456/10854218

如有志同者请发邮箱chenshuai0614@hrbeu.edu.cn联系我！