squre2.m

对于多项式y(i) = b0 + b1 * x + b2 * x^2 + b3 * x^3 + ...... bn * x^n 用递推最小二乘法来估计上式中的b0, b1, ..., bn

% y(i) = b0 + b1 * x + b2 * x^2 + b3 * x^3 + ...... bn * x^n
% 用递推最小二乘法来估计上式中的b0, b1, ..., bn

x = ( -4 : .1 : 4 )';    % 横坐标,以列向量表示
n = size( x, 1 );        % 横坐标的点数
A   = [ ones( n, 1 ), x, x.^2, x.^3, x.^4, x.^5, x.^6, x.^7, x.^8, x.^9, x.^10 ];
% A为如下矩阵:
%    1 x(0) x(0)^2 x(0)^3 x(0)^4 x(0)^5 x(0)^6 x(0)^7 x(0)^8 x(0)^9 x(0)^10
%    1 x(1) x(1)^2 x(1)^3 x(1)^4 x(1)^5 x(1)^6 x(1)^7 x(1)^8 x(1)^9 x(1)^10
%    ......  ......
%    1 x(n) x(n)^2 x(n)^3 x(n)^4 x(n)^5 x(n)^6 x(n)^7 x(n)^8 x(n)^9 x(n)^10

% 生成测试用样本点
y0 = sin( x ); 
y  = y0 + 0.5 * randn( n, 1 );   % 以sin(x)为基础生成一系列随机样本点

Ft = [ 1 1 1 1 1 1 ];
Pt = Ft * Ft';
Kt = Pt * Ft;
bt = zeros( 6, 1 );
I  = eye( 6 );

for i = 1 : 10
    Ft = A( i, 1 : 6 )';
    Kt = Pt * Ft * inv( 1 + Ft' * Pt * Ft);
    bt = bt + Kt * ( y(i) - Ft' * bt );
    Pt = ( I - Kt * Ft') * Pt;
end

figure( 1 ); clf;  
plot( x, y, '.b' ); hold on;    % 显示样本点
yi = A(:, 1:6) * bt;            % 根据估计的系数计算多项式
                                % y(i) = b0 + b1 * x + ... + bn * x^7
disp( sqrt( mean( ( y0 - yi ).^2 ) ) );
disp( bt' );
plot( x, yi, 'g' );             % 绘制 x 为3次的图形