数据拟合线性最小二乘法及其应用(householder变换)

于是 .
最后得到 ，其中 为 阶上三角矩阵， 为 阶矩阵， 为 阶上梯形矩阵.记 ，则 .特别，若 ,
则当m>n时，经n步可将矩阵A化为一个上三角阵，得到 ，
其中R为 阶上三角阵；当m=n时，经n-1步可将A化为一个上三角阵，得到A=QR.
经过以上分析，可以得到矩阵的正交三角化的完整算法.
For j 1 : n
ABSTRACT
This thesis introduces the basic principle of the least square data fitting. According to the least square data fitting ,the paper studies and analyzes the linear least square fitting using the Householder matrix transformations and QR decompositions of the matrix. Using Matlab makes the theory of data fitting easier to be understood. The method of operation and compose of Matlab is quite easy and the result is accurate. It is very convenient for solving practical problems.
此外，上溢和下溢也是计算中需要考虑的问题.当下溢发生时，一些计算机系统自动置其为零，这就可能出现 为零的情形.另外如果x的分量太大，当该分量平方时，便会出现上溢.考虑到对任意的非零实数 单位化向量相同，为了避免溢出现象的出现，我们可用 代替x来构造v（这样做相当于在原来的v之前乘了常数 ）
3.2 QR
x=R(i:n,i);
y=[1;zeros(n-i,1)];
Ht=householder(x,y);
H=blkdiag(eye(i-1),Ht);
Q=Q*H;
R=H*R;
end
5.3
matlab代码如下：
%用QR分解来解方程AX=b
A=[4 2 -1;3 -1 2;11 3 9];%实例演示
b=[2 10 8]';%实例演示
摘
本文介绍了数据拟合中最小二乘法的基本原理，并根据最小二乘数据拟合方法，针对线性最小二乘拟合及其快速算法（运用HOUSEHOLDER变换及QR分解）进行了分析和研究，并通过数学工具Matlab实现了对相关问题的举例应用，使相关数据拟合理论更加生动易懂.
关键词：数据拟合，最小二乘法，householder变换，应用
本文论述了最小二乘数据拟合法的若干原理，介绍了如何对数据应用Matlab软件实现数据拟合，并给出了若干实例帮助理解相关理论及Matlab实现
最小二乘问题在数据拟合、参数估计和函数逼近等方面中有着广泛的应用.
第二章
2.1
给出m组观测数据(ti,yi),i=1,2,…,m.可以认为它是某一模型
得到的，其中 是“真正”的光滑函数， 是误差函数，y(x)是 经 干扰所得到的，因此借用物理学名词，又把 成为白噪音过程， , , ,而 是观测的随机误差.
对于m组观测数据 ，计算广义多项式P(x)在 处所得的值 与 之差的平方和为
目的是要求出 （即确定参数 ）使 达到最小值.由于 是 的函数 .因此问题化为求多元函数 的极小值.由多元函数极值的必要条件
可得n元线代数方程组
其中
n元线代数方程组称为线性最小二乘问题的正规方程组.以 为元素的n元线代数方程组的系数行列式不为零时，可以唯一地解出广义多项式的参数
因Q是正交矩阵，所以 ，若选择b，使得 （2）
那么 将达到极小值，此时由（1）可得 ，且从（1）有 ，故 ，
由（2）可知，n阶上三角形方程组的解b就是最小二乘解，它是非常容易求解的.
这样，可描述正交化方法求最小二乘解的计算步骤为：
1)对 阶矩阵A建立QR分解A=QR；
2)计算 ，形成方程组Rb=c；
3)用回代法解Rb=C得出最小二乘解b．
正交化方法因其数值稳定而显得优于求解法方程，但因对A建立分解的计算量较大，故除非已经有了这种分解的资料可直接引用或者法方程病态，一般还是愿意从法方程来求最小二乘解的.
第五章
5.1
matlab代码如下：
function [H,rho]=householder(x,y)
error('The Column Vectors X and Y Must Have The Same Length!');
end
rho=-sign(x(1))*norm(x)/norm(y);
y=rho*y;
v=(x-y)/norm(x-y);
I=eye(length(x));
H=I-2*v*v';
定理二显示出，对于任意的 都可构造出Householder矩阵H，使Hx的后n-1个分量为零；而且其证明亦告诉我们，可按如下的步骤来构造确定H的单位向量u：
计算 ;
计算
首先，一个自然的问题是，实际计算时， 前的符号如何选取最好.为了使变换后得到的 为负数，则应取
但是这样选取就会出现一个问题，如果x是一个很接近 的向量，计算
第三章
3.1
Gauss消去法是用左乘一系列初等下三角阵来约化一个矩阵为上三角阵.而Householder变换（注释一）是应用非常广泛的另外一种三角化方法.
定义Householder矩阵是指形式为
的矩阵，其中 ，通常记为H.
根据定义可以看出，Householder矩阵H具有良好的性质：
对称性（ ）；
正交性（ ）；
[Q,R]=qr_h(A)
X=R\(Q\b)%AX=b,A=QR
第六章
在房产估价打的线性模型
中， 分别表示税、浴室数目、占地面积、居住面积、车库数目、房屋数目、房龄、建筑类型、户型及壁炉数目，y代表房屋数目，现根据y的28个数据，求出模型中参数的最小二乘结果.
% [Q,R]=qr(A) %调用MATLAB自带的QR分解函数进行验证
% [q,r]=qrhs(A) %调用本函数进行QR分解
% q*r-A %验证A=QR
% q'*q %验证q的正交性
% norm(q) %验证q的标准化，即二范数等于1
%
n=size(A,1);
R=A;
Q=eye(n);
for i=1:n-1
证明见附录.
利用QR分解，就可以实现正交化方法.设 （ ）有线性无关的列， ，并且假定已知A的QR分解.现将Q分块为
并且令
那么
由此即知，x是最小二乘问题的解当且仅当x是 的解.这样一来，最小二乘问题的解可以很容易从上三角方程组 求得.
综合上面的分析，可得求正交化方法的基本步骤为：
（1）计算A的QR分解；
设 ，由于2范数具有正交不变性，故对任意的正交矩阵 有
这样，最小二乘问题
就等价于原最小二乘问题.因此，就可以通过适当选取正交矩阵Q，使原问题转化为较容易求解的最小二乘问题，这就是正交化方法的基本思想.
定理3（QR分解定理）设 ，则A有QR分解：
其中 是具有非负对角元的上三角阵；而且当 且A非奇异时，上述的分解时唯一的.
对称性（ ）；
反射性：对任意的
应用中经常是，确定一个Householder矩阵H的u并不总是单位向量，把它规范化 需要计算 ，因而要求平方根，这是比较费时的.为此，用以下定理的方法：
定理1：设 ，令
则
是一个Householder矩阵.
证明：设 ，那么 ，且
Householder矩阵的一个关键是能把一个给定向量的若干个指定的分量变为零，具体地说就是下列定理：
% 3.如果||x||≠||y||，那么存在HouseHolder矩阵H，使Hx=rho*y，其中rho=-sign(x(1))*||x||/||y||
% 4.Householder矩阵，常用于将一个矩阵A通过正交变换对角阵
%
x=x(:);
y=y(:);
if length(x)~=length(y)
% H*x-rho*y %验证Hx=rho*y
% H'*H %验证正交
%
%关于HouseHolder变换
% 1.H=I-2vv'，其中||v||=1，我们称H为反射(HouseHolder)矩阵，易证H对称且正交
% 2.如果||x||=||y||，那么存在HouseHolder矩阵H=I-2vv，其中v=±(x-y)/||x-y||，使Hx=y
设 是带有n个参数的数据拟合函数，记
决定 中的参数，使
达到极小，这就是最小二乘拟合问题，也常称为最小二乘法或最小二乘逼近，所求得的解称最小二乘解.
2.2
k次代数多项式是 的线性组合
由于 由k+1个系数 唯一确定，所以又称为k+1阶多项式.一般地，设 是x的函数，并且它们是线性无关组，则其线性组合
称为n阶广义多项式，取广义多项式为拟合函数的最小二乘拟合问题称为线性最小二乘问题.
（2）计算 ；
（3）求解上三角方程组 ；
由此可知，实现正交化方法的关键是如何实现矩阵A的QR分解.
3.3
设 ，把 的列向量记作 .
第一步，令 ，利用定理二，取 ，求出
于是
假设进行了 步，得到Householder变换 ，使
其中 是上三角矩阵，假定 ，其中
第k步，令 ，其中 ，其中 ，则 仍为Householder矩阵，
定理2：设 ， ，且假定 ，则可找到一个Householder矩阵H使
其中
证明：作 ，求出 ，则 即为所求的Householder矩阵.
因为
所以
上述定理的证明是构造性的，即它叙述了计算u与 的过程，一旦 的正负号确定后，就可以求出
以及
如果 异号，在计算 时就会发生抵消，因此，我们取 .
这样我们可以运用Householder变换把系数矩阵化为上梯形矩阵.
时，就会出现两个相近的数相减，而导致严重地损失有效数字，这里 分别表示向量v和x的第一个分量.不过，幸运的是，只要对上式做一简单的等价变形，就可避免这一问题的出现，事实上，注意到
只要在 时使用上面式子Βιβλιοθήκη Baidu计算 ，就会避免出现两个相近的数相减的情形.
其次，注意到
其中 ，我们就没有必要非求出u不可，而只需求出 即可.然而在实际计算时，将 为第一个分量为1的向量是方便的，这是因为这样正好可以把 个分量保存在 个化为0的分量位置上，而 的第一个分量1就无需保存了.
特别地，当 时，相应的n元线性代数方程组为
其中
2.3 小结
最小二乘问题的解归结为求正则方程组.由于正则方程系数矩阵的条件数是原矛盾方程组系数条件数的平方，正则方程的病态程度大大增加.因此，有必要采用有较高数值稳定性的计算方法.一般使用H-变换，平方根法和Gram-Schmidt正交化法.这三种方法用于一次性处理时，H-变换所需的运算量较少，而得到最广泛的使用.平方根法便于形成递推算法而首先被用于最小二乘估计的递推算法中，但与直接求解正则方程的递推算法相比较，运算量有所增加.
5.2
matlab代码如下：
function [Q,R]=qr_h(A)
%基于Householder变换，将方阵A分解为A=QR，其中Q为正交矩阵，R为上三角阵
%
%参数说明
% A：需要进行QR分解的方阵
% Q：分解得到的正交矩阵
% R：分解得到的上三角阵
%
%实例说明
% A=[-12 3 3;3 1 -2;3 -2 7];
Key Words:data fittingleast square methodhouseholder transformation application
第一章
实验数据的正确处理关系到能否达到实验目的、得出明确结论.在进行实验数据的处理分析时，数据拟合是经常采用的方法之一.数据拟合的目标是寻找一条光滑数据，使之在某种准则下最佳地拟合数据.
If j m
end
end
function :
n length(x)
else
else
end
end
第四章
求解线性最小二乘问题的快速算法，主要是避免构造法方程组，而直接从矛盾方程组Ab=Y入手．常应用Householder变换把系数矩阵A正交三角化，使 ，其中R为n阶上三角阵，0为 的零矩阵，Q是一个m阶正交矩阵．并把m维向量QY相应地分块n维向量c与 维向量d，即 .于是 （1）
%求解正交对称的Householder矩阵H，使Hx=rho*y，其中rho=-sign(x(1))*||x||/||y||
%
%参数说明
% x：列向量
% y：列向量，x和y必须具有相同的维数
%
%实例说明
% x=[3 5 6 8]';
% y=[1 0 0 0]';
% [H,rho]=householder(x,y);