用正交多项式实现最小二乘拟合

1. 目的
   用最小二乘法求一个多项式，基于残差平方和最小原则拟合下列数据点
   

熟悉最小二乘拟合方法，掌握正交多项式的构造方法，利用正交多项式实现最小二乘拟合，可以避免求解正规方程组。实现MATLAB程序设计，求出拟合曲线，均方误差，并画出拟合函数的图像。

1. 理论
   
2. 程序

最小二乘拟合函数：
function y = least_squares_fitting_fun(xi,yi,wi,n)
x = sym('x');     %打印拟合函数中用到的字符'x'
m = length(xi);
P=zeros(n,m);     %构建矩阵用于存放原始点处P的值
 for i=1:m        %首先将P的第一行赋值为1,即xi(i)处的p0值都为1
    P(1,i)=1;
end
for k=2:n+1       %根据公式计算α，β
    a(k)=sum(wi(i).*xi.*P(k-1,:).^2)/sum(wi(i).*P(k-1,:).^2);
    if k>2
        b(k)=sum(wi(i).*P(k-1,:).^2)/sum(wi(i).*P(k-2,:).^2);
    end
    for i=1:m     %计算在数据点xi(i)处Pk处的值
        if k>2
        P(k,i)=(xi(i)-a(k))*P(k-1,i)-b(k)*P(k-2,i);
        else
           P(k,i)=(xi(i)-a(k))*P(k-1,i);
        end
    end
end
for k=1:n+1       %计算a*
    ak(k)=sum(wi(i).*yi.*P(k,:))/sum(wi(i).*P(k,:).^2);
end
p=cell(n,1);       %利用元胞数组存放P（里面是符号表达式）
p{1}=1;          %将P0赋值为1
for k=2:n+1       %根据α，β和a求出P1到Pk的表达式
    if k>2
        p{k}=(x-a(k))*p{k-1}-b(k)*p{k-2};
    else
        p{k}=(x-a(k))*p{k-1};
    end
end
y = 0;
for k=1:n+1       %得出最终的拟合函数表达式
    y=y+ak(k)*p{k};
end
                 %计算平方误差d
Aak = 0.0;
for k = 1:n+1
    Aak = Aak+sum(P(k,:).^2)*ak(k)^2;
end
d = sum(yi.^2)-Aak;
sprintf('均方误差为：%f',d)
                  %画图
plot(xi,yi,'*'); hold on;
first_a = min(xi);
end_b = max(xi);
fplot(y, [first_a-1,end_b+1]);
title('最小二乘拟合函数图像');
end
main函数：
xi = [-1.0 -0.5 0.0 0.5 1.0 1.5 2.0];       
yi = [-4.447 -0.452 0.551 0.048 -0.447 0.549 4.552];
wi = [1 1 1 1 1 1 1];         %数据的权重
n = 3;                    %拟合函数阶数
least_squares_fitting_fun(xi,yi,wi,n);

4. 实验结果