系统辨识最小二乘法大作业
哈工大系统辨识实验三
zl=[z(3:16)]'
hl1=-z(2:15)'
hl2=-z(1:14)'
hl3=u(2:15)'
hl4=u(1:14)'
Hl=[hl1,hl2,hl3,hl4];
theta=((Hl'*Hl)^-1)*Hl'*zl;
disp('辨识后得到的参数：')
disp(theta)
end
function m=m_xulie(n,k)%n为M序列的阶数k为序列的长度
………………………………………………………………………………………….....（2）
设
，
则式（2）可写为
…………………………………………………………………………….（3）
式中：y为N维输出向量； 为N维噪声向量； 为 维参数向量； 为 测量矩阵。为了尽量减小噪声 对 估值的影响，应取 ，即方程数目大于未知数数目。
通过本次上机实验课程，我体会到了最小二乘法辨识是一种近代的系统辨识方法，其优点是只需要考虑系统的输入输出特性，而不强调系统的内部机理。适用于单输入单输出系统的差分方程模型中的参数估计。
最后感谢老师在指导我们上机实验当中所付出的努力。
当系统中没有白噪声时，理论上应该得到准确的参数，我们把模型中的v(k)项剔除，得到以下的辨识结果：
这说明当不存在白噪声时，可以精确获得精确的估计值 ，估计值 的均方差随着噪声的均方差的增大而增大。
8
在这次上机实验当中，我们跟着老师学会了运用最小二乘法实现对系统的参数的估计，估计结果与系统的实际参数存在一定的偏差，原因是白噪声对系统的干扰造成的。
x(j)=X(n);
end
x=-2*x+ones(1,k);
系统辨识相关分析最小二乘
相关分析法辨识系统单位脉冲响应1辨识原理对于下图示的单输入单输出线性系统，其输入输出的因果关系可用卷积公式描述。
公式为：0()()()y t g x t d λλλ∞=-⎰把变量t 换成t +τ，上式两边同乘以x (t )，取时间的平均值，得11lim()(+)(){lim()(+)}22TTTTT T x t y t dt g x t x t dt d TTτλτλλ∞--→∞→∞=-⎰⎰⎰即 0()()()x y x R g R d τστλλ∞=-⎰上式即为维纳-霍夫方程，其给出了输入的自相关函数，输入、输出的互相关函数及脉冲响应函数三者之间的关系。
令x (t )为白噪声信号，则其自相关函数为：()(), ()()x x R k R k τδττλδτλ=-=-代入维纳-霍夫方程得：()()()()xy x R g R d kg τλτλλτ∞=-=⎰则有：()()xy R g kττ=这样，只要记录x(t)、y(t)的值，并计算它们的互相关函数，即可求得脉冲响应函数g(τ)。
在系统有正常输入的情形下，辨识脉冲响应的原理图如下图所示。
2辨识过程2.1预备实验以二阶系统22()2G s s s ++=作为辨识对象。
在实验前首先要进行预备实验，以了解系统特性。
通过简单阶跃响应确定系统过度过程时间T s 大约为11s ，如下图所示。
给系统施加不同周期的正弦信号，系统输出为输入的0.707倍时，确定截止频率f M 大约为0.318Hz 。
2.2选择二位式伪随机序列的参数（1）选择t 和N 由0.3Mt f ∆≤，得0.94t s ∆≤。
因为系统的时间常数1nT s ζω=，根据时间常数可按照0.050.1t T ∆= （）选择t ∆。
由二位式伪随机序列周期要大于系统过渡过程时间，若t ∆选择0.94s ，则由(1)s N t T -⨯∆≥，得12.7021N ≥；若t ∆选择0.195s ，则由(1)s N t T -⨯∆≥，得57.4103N ≥。
系统辨识作业
1、利用cftool拟合曲线得Coefficients (with 95% confidence bounds):a = 0.4601 (0.3711, 0.5492)b = -5.1 (-6.108, -4.092)c = 13.42 (11, 15.83)2、最小二乘法：最小二乘的思想就是寻找一个θ 的估计值θ ，使得各次测量的Z i (i = 1,⋯ m)与由估计θ 确定的量测估计Z i = H iθ 之差的平方和最小，由于此方法兼顾了所有方程的近似程度，使整体误差达到最小，因而对抑制误差是有利的。
但是最小二乘一般方法的估计精度不够高，这是由于对各个测量数据同等对待，而各次测量数据一般不会在相同的条件下获得，造成测量数据的置信度不变较大，当新数据源源而来时，将出现以下问题：1.数据增加，要求计算机的存储空间增加2、每增加一组数据，即作一次求逆，导致计算量增加，难以用于在线辨识。
只有当模型的噪声项是独立的随机变量时，普通最小二乘法才能得到真实参数的无偏估计值，否则所得到的估计值是有偏的。
因此，最小二乘法有很多改进算法，虽然没有一个是完美的，但是能够适应不同的情况、条件，对应选择不同的算法，其各自的性能及优缺点如下：广义最小二乘法的优点是计算精度高，估计的效果比较好，是无偏估计，但广义最小二乘法缺点是计算量大，其收敛是比较缓慢的，为了得到准确的参数估值，往往需要进行多次迭代计算，另外，对于循环的循环性没有给出证明，并非总是收敛于最优估值上。
由于一般情况下，系统信噪比比较低，准则函数为非单值函数（即存在多个局部极小值），如果初值给的不合理，用GLS方法得到的将是局部极小值，若想得到总体最优解，初值应接近该最优值。
递推算法实现了实时控制，减少了计算量和存储量，但未解决最小二乘法的递推算法有偏估计问题。
矩形窗（限定记忆）RLS 方法需要保留一定的数据存储量，此存储量大小取决于矩形窗宽度，因而在应用范围上有一定程度的限制。
系统辨识之最小二乘法
系统辨识之最小二乘法方法一、最小二乘一次性算法：首先对最小二乘法的一次性辨识算法做简要介绍如下：过程的黑箱模型如图所示：其中u(k)和z(k)分别是过程的输入输出，)(1-z G 描述输入输出关系的模型，成为过程模型。
过程的输入输出关系可以描述成以下最小二乘格式：)()()(k n k h k z T +=θ （1）其中z(k)为系统输出，θ是待辨识的参数，h(k)是观测数据向量，n(k)是均值为0的随机噪声。
利用数据序列{z （k ）}和{h （k ）}极小化下列准则函数：∑=-=Lk T k h k z J 12])()([)(θθ （2）使J 最小的θ的估计值^θ，成为最小二乘估计值。
具体的对于时不变SISO 动态过程的数学模型为 )()()()()(11k n k u z B k z z A +=-- （3）应该利用过程的输入、输出数据确定)(1-z A 和)(1-Z B 的系数。
对于求解θ的估计值^θ，一般对模型的阶次a n ，b n 已定，且b a n n >；其次将（3）模型写成最小二乘格式)()()(k n k h k z T +=θ （4）式中=------=T n n T b a b a b b b a a a n k u k u n k z k z k h ],,,,,,,[)](,),1(),(,),1([)(2121 θ （5）L k ,,2,1 =因此结合式（4）（5）可以得到一个线性方程组L L L n H Z +=θ （6）其中==T L TL L n n n n L z z z z )](),2(),1([)](),2(),1([ （7）对此可以分析得出，L H 矩阵的行数为),max(b a n n L -，列数b a n n +。
在过程的输入为2n 阶次，噪声为方差为1，均值为0的随机序列，数据长度)(b a n n L +>的情况下，取加权矩阵L Λ为正定的单位矩阵I ，可以得出：L T L L T L z H H H 1^)(-=θ （8）其次，利用在Matlab 中编写M 文件，实现上述算法。
机械系统辨识及仿真最小二乘法
作业四宋家亮15030024一、实验题目用MatLAB辨识系统，系统输出分别叠加两种不同类型的噪声，针对每一种情况使用最小二乘整批算法、递推算法和广义最小二乘法实现辨识。
图1 图2二、实验思路及目的对图1、图2分别通过构造的系统给定输入（白噪声和M序列）并叠加噪声（白噪声和有色噪声），测得输出数据并显示，对未知系统定阶并用整批、递推、广义最小二乘算法进行辨识，显示辨识结果；比较三种算法对不同噪声模型的辨识精度（λ取相同的值），显示辨识结果的脉冲响应图像并于理想系统响应对比，对结果给予合理的解释；改变 值(0.01,0.02,0.03等)，比较辨识结果的精度，说明信噪比对各种算法辨识精度的影响。
三、实验过程1、设计界面利用guide编辑器针对所要实现的功能添加控件进行界面初步的设计，设计结果如图1.1所示图1.12、设置控件属性对界面的控件设置属性，形成实验所需要的最终界面，如图2.1所示：图2.1其中系统噪声和输入信号设置成可选择的属性，其value值作为后面函数中的控制条件，这个刚开始走了弯路，直接用string属性做控制条件，结果发现根本实现不了，请教大神后选择value值作为控制条件，因为value值默认为1，当通过下拉选项选择下一个信号时其value值会自动加1。
如图2.2所示。
图2.23、设计guide回调函数（1）噪声参数输入的回调函数本实验设计成噪声参数可输入的程序，所以为达到这一目的，添加了参数输入控件，通过如图3.1所示打开回调函数，添加如下程序，可以将输入的容转变成数字用于后续函数的运算。
global zaoshengxishuzaoshengxishu=get(handles.zaoshengxishu,'String'); zaoshengxishu=str2num(zaoshengxishu);%将输入的字符串λ变成数值图3.1（2）开始辨识按钮的回调函数（主体程序）如图3.2打开开始辨识按钮回调函数，添加主程序（见附录），这样通过点击该按钮可实现本实验的所有功能。
系统辨识最小二乘法大作业 (2)
系统辨识大作业最小二乘法及其相关估值方法应用学院：自动化学院学号：姓名：日期：基于最小二乘法的多种系统辨识方法研究一、实验原理1.最小二乘法在系统辨识中用得最广泛的估计方法是最小二乘法(LS)。
设单输入-单输出线性定长系统的差分方程为(5.1.1)式中:为随机干扰；为理论上的输出值。
只有通过观测才能得到，在观测过程中往往附加有随机干扰。
的观测值可表示为(5.1.2)式中:为随机干扰。
由式(5.1.2)得(5.1.3)将式(5.1.3)带入式(5.1.1)得(5.1.4)我们可能不知道的统计特性，在这种情况下，往往把看做均值为0的白噪声。
设(5.1.5)则式(5.1.4)可写成(5.1.6)在观测时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。
因此假定不仅包含了的测量误差，而且包含了的测量误差和系统内部噪声。
假定是不相关随机序列(实际上是相关随机序列)。
现分别测出个随机输入值，则可写成个方程，即上述个方程可写成向量-矩阵形式(5.1.7) 设则式(5.1.7)可写为(5.1.8)式中：为维输出向量；为维噪声向量；为维参数向量；为测量矩阵。
因此式(5.1.8)是一个含有个未知参数，由个方程组成的联立方程组。
如果，方程数少于未知数数目，则方程组的解是不定的，不能唯一地确定参数向量。
如果，方程组正好与未知数数目相等，当噪声时，就能准确地解出(5.1.9)如果噪声，则(5.1.10)从上式可以看出噪声对参数估计是有影响的，为了尽量较小噪声对估值的影响。
在给定输出向量和测量矩阵的条件下求系统参数的估值，这就是系统辨识问题。
可用最小二乘法来求的估值，以下讨论最小二乘法估计。
2.最小二乘法估计算法设表示的最优估值，表示的最优估值，则有(5.1.11)写出式(5.1.11)的某一行，则有(5.1.12) 设表示与之差，即-(5.1.13)式中成为残差。
把分别代入式(5.1.13)可得残差。
设则有(5.1.14) 最小二乘估计要求残差的平方和为最小，即按照指数函数(5.1.15) 为最小来确定估值。
安大系统辨识作业
2012系统辨识实验姓名：周自飞学号：P4*******年级专业：09自动化线性系统参数估计的最小二乘法一、实验目的：1、掌握线性离散系统的数学模型；2、掌握线性离散系统在无噪声时的单位脉冲响应；3、掌握线性系统参数估计的最小二乘法。
二、实验原理：1、基本最小二乘算法2、递推最小二乘算法3、渐消记忆递推最小二乘算法4、辅助变量递推最小二乘算法5、增广递推最小二乘算法三、实验内容：设带有噪声的离散系统模型为：()()()()()()k213.0-13.0-+.028.11-412+ukkukyk--yξky=现给定系统的一列长度为50的输入序列：u(k)= 1,0.3,-0.5,0.9,-0.5,-0.3,-0.2,0.4,0.8,-0.6,-0.1,0,0.1,0.5,-0.6,-0.2,0.3,0.9,0.5,0.2,-0.6,-0.3,-0.1,0.7,0,0.3,-0.6,0.3,0.1,0.5,-0.7,-0.4,-0.9,-0.6,-0.2,-0.4,-0.2,0.1,-0.1,0.1,0.9,0.5,0.3,0.7,0.4,-0.2,-0.7,-0.2,0.1,01. 构造离散系统模型函数，给出无噪声时系统在单位脉冲输入下的响应序列；%zhouzifei.mb=[0.3 -0.213];a=[1 -1.28 0.41];impz(b,a,30);title(‘系统单位脉冲响应’);axis([-1 30 -0.2 0.5]);2. 给出系统在给定输入u(k)和噪声()kξ下的输出信号y(k)；%zifei2.m3. 利用已有输入信号u(k)和在第2步中得出的输出信号y(k)，使用下列方法辨识系统的参数：(1) 基本最小二乘算法%LS.m参数：a1 =-1.2803a2 =0.4138b0 = 0.2867b1 = -0.2104(2) 递推最小二乘算法%RLS.m参数：a1 =-1.2800 a2 =0.4100 b0 = 0.3000 b1 = -0.2130(3) 渐消记忆递推最小二乘算法%FMRLS参数：a1 =-1.2800 a2 =0.4100 b0 = 0.3000 b1 = -0.2130(4) 辅助变量递推最小二乘算法%IVLS.m参数：a1 =-1.2800 a2 =0.4100b0 = 0.3000 b1 = -0.2130(5) 增广递推最小二乘算法（可选做，有加分）%ARLS.m参数：a1 =-1.2800 a2 =0.4100 b0 = 0.3000 b1 = -0.2130四、实验结果：最小二乘法思想是使各次观测值和计算值之间差值的平方乘以度量其精确度的数值后的和为最小。
系统辨识大作业1201张青
《系统辨识》大作业学号：********班级：自动化1班姓名：**信息与控制工程学院自动化系2015-07-11第一题模仿index2，搭建对象，由相关分析法，获得脉冲响应序列ˆ()g k，由ˆ()g k，参照讲义，获得系统的脉冲传递函数()G z和传递函数()G s;应用最小二乘辨识,获得脉冲响应序列ˆ()g k;同图显示两种方法的辨识效果图；应用相关最小二乘法，拟合对象的差分方程模型；构建时变对象，用最小二乘法和带遗忘因子的最小二乘法，(可以用辨识工具箱) 辨识模型的参数，比较两种方法的辨识效果差异;答：根据index2搭建结构框图：相关分析法：利用结构框图得到UY 和tout其中的U就是题目中要求得出的M序列，根据结构框图可知序列的周期是1512124=-=-=npN。
在command window中输入下列指令，既可以得到脉冲相应序列()g k：aa=5;NNPP=15;ts=2; RR=ones(15)+eye(15); for i=15:-1:1UU(16-i,:)=UY(16+i:30+i,1)'; endYY=[UY(31:45,2)];GG=RR*UU*YY/[aa*aa*(NNPP+1)*ts]; plot(0:2:29,GG) hold onstem(0:2:29,GG,'filled') Grid;title('脉冲序列g(τ)')最小二乘法建模的响应序列由于是二阶水箱系统，可以假设系统的传递函数为221101)(sa s a sb b s G +++=，已知)(τg ，求2110,,,a a b b已知G （s ）的结构，用长除法求得G(s)的s 展开式，其系数等于从 )( g 求得的各阶矩，然后求G(s)的参数。
得到结果： a1 =-1.1561 a2 =0.4283 b0 =-0.0028 b1=0.2961在command window 中输入下列指令得到传递函数：最小二乘一次算法相关参数%最小二乘法一次完成算法 M=UY(:,1); z=UY(:,2); H=zeros(100,4); for i=1:100 H(i,1)=-z(i+1); H(i,2)=-z(i); H(i,3)=M(i+1); H(i,4)=M(i); endEstimate=inv(H'*H)*H'*(z(3:102)) %结束得到相关系数为：Estimate =-0.7866 0.1388 0.5707 0.3115带遗忘因子最小二乘法：%带遗忘因子最小二乘法程序M=UY(:,1);z=UY(:,2);P=1000*eye(5); %Theta=zeros(5,200); %Theta(:,1)=[0;0;0;0;0];K=zeros(4,400); %K=[10;10;10;10;10];lamda=0.99;%遗忘因数for i=3:201h=[-z(i-1);-z(i-2);M(i);M(i-1);M(i-2)];K=P*h*inv(h'*P*h+lamda);Theta(:,i-1)=Theta(:,i-2)+K*(z(i)-h'*Theta(:,i-2));P=(eye(5)-K*h')*P/lamda;endi=1:200;figure(1)plot(i,Theta(1,:),i,Theta(2,:),i,Theta(3,:),i,Theta(4,:),i,Theta(5,:) )title('带遗忘因子最小二乘法')grid%结束Estimate 可由仿真图得出，可知两种方法参数确定十分接近。
系统辨识最小二乘法大作业
系统辨识大作业最小二乘法及其相关估值方法应用学号：2012302259姓名：王家琦基于最小二乘法的多种系统辨识方法研究1.最小二乘法的引出在系统辨识中用得最广泛的估计方法是最小二乘法(LS)。
设单输入-单输出线性定长系统的差分方程为，，，，(5.1.1)式中:为随机干扰；为理论上的输出值。
只有通过观测才能得到，在观测过程中往往附加有随机干扰。
的观测值可表示为(5.1.2)式中:为随机干扰。
由式(5.1.2)得(5.1.3)将式(5.1.3)带入式(5.1.1)得(5.1.4)我们可能不知道的统计特性，在这种情况下，往往把看做均值为0的白噪声。
设(5.1.5)则式(5.1.4)可写成(5.1.6)在观测时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。
因此假定不仅包含了的测量误差，而且包含了的测量误差和系统内部噪声。
假定是不相关随机序列(实际上是相关随机序列)。
现分别测出个随机输入值，，，，，，，，则可写成个方程，即上述个方程可写成向量-矩阵形式(5.1.7)设，，则式(5.1.7)可写为(5.1.8)式中：为维输出向量；为维噪声向量；为维参数向量；为测量矩阵。
因此式(5.1.8)是一个含有个未知参数，由个方程组成的联立方程组。
如果，方程数少于未知数数目，则方程组的解是不定的，不能唯一地确定参数向量。
如果，方程组正好与未知数数目相等，当噪声时，就能准确地解出(5.1.9)如果噪声，则(5.1.10) 从上式可以看出噪声对参数估计是有影响的，为了尽量较小噪声对估值的影响。
在给定输出向量和测量矩阵的条件下求系统参数的估值，这就是系统辨识问题。
可用最小二乘法来求的估值，以下讨论最小二乘法估计。
2.最小二乘法估计算法设表示的最优估值，表示的最优估值，则有(5.1.11)，写出式(5.1.11)的某一行，则有，，，，(5.1.12)设表示与之差，即-，，，，(5.1.13)式中成为残差。
把，，，分别代入式(5.1.13)可得残差，，，。
系统辨识与参数估计大作业
系统辨识与参数估计大作业第一题 递推最小二乘估计参数考虑如上图所示的仿真对象，选择模型结构为：)()2()1()2()1()(2121k v k u b k u b k z a k z a k z +-+-=-+-+，其中)(k v 是服从)1,0(N 正态分布的不相关随机噪声；输入信号)(k u 采用4阶逆M 序列，特征多项式取41)(s s s F ⊕⊕=，幅度为1，循环周期为bit N p 62=；控制λ值，使数据的噪信比分别为10％，73％，100％三种情况。
加权因子1)(=Λk ；数据长度L=500;初始条件取001.0)0(ˆ,10)0(6==θI P ， （1） 利用递推最小二乘算法在线估计参数，（2） 利用模型阶次辨识方法（AIC 准则），确定模型的阶次。
（3） 估计噪声)(k v 的方差和模型静态增益K （4） 作出参数估计值随时间的变化图 答：设过程的输入输出关系可以描述成()()()Tz k h k n k θ=+()z k 是输出量，()h k 是可观测的数据向量，n (k)是均值为0的随机噪声[]()(1),(2),(1),(2)Th k z k z k u k u k =------[]1212,,,Ta ab b θ=选取的模型为结构是1212()(1)(2)(1)(2)z k a z k a z k bu k bu k =----+-+- 12121.5,0.7, 1.0,0.5a a b b =-===加权最小二乘参数估计递推算法RWLS 的公式如下，11()(1)()()(1)()()()(1)()()()(1)()()()(1)T T TK k p k h k h k p k h k k k k K k z k h k k p k I K k h k p k θθθ-⎡⎤=--+⎢⎥Λ⎣⎦⎡⎤=-+--⎣⎦⎡⎤=--⎣⎦为了把p(k)的对称性，可以把p(k)写成1()(1)()()()(1)()()T T p k p k K k K k h k p k h k k ⎡⎤=---+⎢⎥Λ⎣⎦如果把()k Λ设成1的时候，加权最小二乘法就退化成最小二乘法。
系统辩识作业题
系统辨识大作业
一.设SlSO系统差分方程为
y(k)=—α1y(k-1)-a2y(k-2)+bλu(k-1)+b2u(k-2)+ξ{k)
辨识参数向量为θ=[q a2b l b2]r,输入输出数据详见数据文件UyLtXt—uy3.txtoξ(k)为噪声方差各异的白噪声或有色噪声。
试求解：
1)用n元一次方程解析法，再求其平均值方法估计。
2)用最小二乘及递推最小二乘法估计。
；
3)用辅助变量法及其递推算法估计
4)用广义最小二乘法及其递推算法估计
5)用夏氏偏差修正法、夏氏改良法及其递推算法估计
6)用增广矩阵法估计
7)分析噪声父攵)特性；
二.用极大似然法估计6。
三.以上题的结果为例，进行：
1.分析比较各种方法估计的精度；
2.分析其计算量；
3.分析噪声方差的影响；
4.比较白噪声和有色噪声对辨识的影响。
四.系统模型阶次的辨识：
1.用三种方法确定系统的阶次并辨识；
2.分析噪声对定阶的影响；
3.比较所用三种方法的优劣及有效性；
五.给出由正弦输入求取系统开环频率响应特性曲线的辨识方法。
六.提出一种自己创造的辨识新方法，并用所给数据进行辨识验证。
注：闭卷考试时提交大作业报告。
系统辨识作业_梯度法最小二乘法
系统辨识一、最小二乘算法仿真对象为：)()2(5.0)1()2(7.0)1(5.1)(k v k k u k k z k z +-+-=-+-+，)(k v 为服从N(0,1)分布的白噪声，输入信号采用4阶M 序列。
辨识模型为：)()2()1()2()1()(2121k v k z b k z b k z a k z a k z +-+-+----=，观测数据长度取400=L ，加权矩阵取I =Λ。
则matlab 程序如下:% M array constructionNp=15;r=4;X1=1;X2=1;X3=1;X4=1;m_length = r*Np;a=1;for i=1:1:m_lengthY4=X4;Y3=X3;Y2=X2;Y1=X1;X4=Y3;X3=Y2;X2=Y1;X1=xor(Y3,Y4);if Y4==0M(i)=-a;elseM(i)=a;endendfigure;i=1:1:m_length;plot(i,M);% 白噪声noise = zeros(1,m_length);for i=1:1:m_lengthtemp = noise + 0.5*rands(1,m_length);noise = temp;endnoise = noise/12;%noise = temp;% parameter of systemn=2;d=1;a1=-1;a2=0.5;b1=1;b2=0.5;S_U0=0.2;S_Y0=0.2;% generate u,yu0=ones(1,m_length)*S_U0;U = M + u0 + noise;figure;i=1:1:m_length;plot(i,U);%formulationy(1)=0;y(2)=0;y(3)=0;Y(1)=S_Y0+y(1)+noise(1);Y(2)=S_Y0+y(2)+noise(2);Y(3)=S_Y0+y(3) +noise(3); for k=4:m_lengthy(k) = b1*U(k-1-d)+b2*U(k-2-d)-a1*y(k-1)-a2*y(k-2);Y(k)=S_Y0+y(k)+noise(k);endfigure;i=1:1:m_length;plot(i,Y);%辨识% premitive valuec=2;P = (c^3)*eye(3*n+d);sita(:,3) = [0,0,0,0,0,0,0]';alf = 0.995;% M%compute U0,Y0sum_U = 0;sum_Y = 0;for k=1:1:m_lengthsum_U = sum_U + U(k);sum_Y = sum_Y + Y(k);endt_U0 = sum_U/m_length;t_Y0 = sum_Y/m_length;for k=1:1:m_lengtht_u(k) = U(k) - t_U0;t_y(k) = Y(k) - t_Y0;endfigure;i=1:1:m_length;plot(i,t_u);figure;i=1:1:m_length;plot(i,t_y);v(1)=0;v(2)=0;v(3)=0;for k=4:1:m_lengthfai = [-t_y(k-1),-t_y(k-2),t_u(k-1-d),t_u(k-2-d),v(k-1),v(k-2),v(k-3)]';v(k) = t_y(k) - fai'*sita(:,k-1);W = P*fai;dan = 1/(alf + fai'*W);sita(:,k) = sita(:,k-1) + dan*W*v(k);P = ( P - dan * W*W')/alf;end仿真结果如下图所示：则各个参数的估计值为：⎪⎪⎪⎪⎪⎭⎫⎝⎛-=⎪⎪⎪⎪⎪⎪⎭⎫⎝⎛∧∧∧∧5048.00062.17003.04961.1121bbaa二、梯度校正辨识方法仿真对象为)()2(5.0)1()2(7.0)1(5.1)(kvkkukkzkz+-+-=-+--，其中)(kv为服从N(0,1)分布的白噪声，取初值0)0(=∧θ，1=α，1.0=c。
系统辨识最小二乘法大作业
系统辨识大作业最小二乘法及其相关估值方法应用学号：2012302259姓名：王家琦基于最小二乘法的多种系统辨识方法研究1.最小二乘法的引出在系统辨识中用得最广泛的估计方法是最小二乘法(LS)。
设单输入-单输出线性定长系统的差分方程为，，，，(5.1.1)式中:为随机干扰；为理论上的输出值。
只有通过观测才能得到，在观测过程中往往附加有随机干扰。
的观测值可表示为(5.1.2)式中:为随机干扰。
由式(5.1.2)得(5.1.3)将式(5.1.3)带入式(5.1.1)得(5.1.4)我们可能不知道的统计特性，在这种情况下，往往把看做均值为0的白噪声。
设(5.1.5)则式(5.1.4)可写成(5.1.6)在观测时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。
因此假定不仅包含了的测量误差，而且包含了的测量误差和系统内部噪声。
假定是不相关随机序列(实际上是相关随机序列)。
现分别测出个随机输入值，，，，，，，，则可写成个方程，即上述个方程可写成向量-矩阵形式(5.1.7)设，，则式(5.1.7)可写为(5.1.8)式中：为维输出向量；为维噪声向量；为维参数向量；为测量矩阵。
因此式(5.1.8)是一个含有个未知参数，由个方程组成的联立方程组。
如果，方程数少于未知数数目，则方程组的解是不定的，不能唯一地确定参数向量。
如果，方程组正好与未知数数目相等，当噪声时，就能准确地解出(5.1.9)如果噪声，则(5.1.10) 从上式可以看出噪声对参数估计是有影响的，为了尽量较小噪声对估值的影响。
在给定输出向量和测量矩阵的条件下求系统参数的估值，这就是系统辨识问题。
可用最小二乘法来求的估值，以下讨论最小二乘法估计。
2.最小二乘法估计算法设表示的最优估值，表示的最优估值，则有(5.1.11)，写出式(5.1.11)的某一行，则有，，，，(5.1.12)设表示与之差，即-，，，，(5.1.13)式中成为残差。
把，，，分别代入式(5.1.13)可得残差，，，。
系统辨识大作业
系统辨识大作业专业班级：自动化09-3学号：09051325姓名：吴恩作业一:设某物理量Y与X满足关系式2=++，实验获得一批数据Y aX bX c如下表，试辨识模型参数,,a b c。
（15分）解答：问题描述：由题意知，这是一个已知模型为Y=aX2+bX+c，给出了10组实验输入输出数据，要求对模型参数a，b，c进行辨识。
问题求解：这里对该模型参数辨识采用最小二乘法的一次算法（LS）求解。
2=++可以写成矩阵形式Y=AE+e;其中A=[X^2,X,1]构成, Y aX bX c利用matlab不难求解出结果。
运行结果：利用所求的的参数，求出给定的X对应的YE值，列表如下做出上表的图形如下12345678910xyy=ax 2+bx+c 参数求解结果分析：根据运行结果可以看出，拟合的曲线与真是观测的数据有误差，有出入，但是误差较小，可以接受。
出现误差的原因，一方面是由于给出的数据只有十个点，数据量太少，难以真正的充分的计算出其参数，另外，该问题求解采用的是LS 一次算法，因此计算方法本身也会造成相应的误差。
作业二：模仿实验二，搭建对象，由相关分析法，获得脉冲相应序列()g k，由()G z；和传递函数g k，参照讲义，获得系统的脉冲传递函数()G s及应用相关最小二乘法，拟合对象的差分方程模型；加阶跃()扰动，用最小二乘法和带遗忘因子的最小二乘法，辨识二阶差分方程的参数，比较两种方法的辨识差异；采用不少于两种定阶方法，确定对象的阶次。
对象模型如图：利用相关分析法，得到对象的脉冲相应序列。
如下图：（1）.由脉冲相应序列，求解系统的脉冲传递函数G（z）Transfer function:0.006072 z^2 + 0.288 z + 0.1671-------------------------------z^2 + 0.1018 z - 0.7509Sampling time: 2(2).由脉冲相应序列求解系统的传递函数G(s)Transfer function:(0.04849+2.494e-018i)-----------------------s^2 + 0.1315 s + 0.6048(3).利用相关最小二乘法拟合系统的差分方程模型如下：(4).在t=100，加入一个0.5的阶跃扰动，，利用RLS求解差分方程模型：RLS加入遗忘因子之后与未加之前的曲线情况如下：未加遗忘因子之前参数以及残差的计算过程加入0.99的遗忘因子得到的参数辨识过程与残差的变化过程根据上面两种方法所得到的误差曲线和参数过渡过程曲线，我们可以看出来利用最小二乘法得到的参数最终趋于稳定，为利用带遗忘因子的最小二乘算法，曲线参数最终还是有小幅度震荡。
系统辨识大作业.
一、 问题描述考虑仿真对象：()0.9(1)0.15(2)0.02(3)0.7(1)0.15(2)()z k z k z k z k u k u k e k +-+-+-=---+ e() 1.0e(1)0.41e(2)(),~(0,1)k k k v k v N λ+-+-=式中，u(k)和z(k)是输入输出数据，v(k)是零均值、方差为1的不相关的随机噪声；u(k)采用与e(k)不相关的随机序列。
1. 设计实验，产生输入输出数据；2. 使用基本最小二乘法估计参数；3. 考虑其他适用于有色噪声的辨识方法估计参数；4. 模型验证。
二、 问题分析对于单输入单输出系统（Single Input Single Output, SISO ），如图 1所示，将待辨识的系统看作“灰箱”，它只考虑系统的输入输出特性，而不强调系统的内部机理。
图 1中，输入u(k)和输出z(k)是可以测量的，1()G z -是系统模型，用来描述系统的输入输出特性，y(k)是系统的实际输出。
1()N z -是噪声模型，v(k)是均值为零的不相关随机噪声，e(k)是有色噪声。
图 1 SISO 系统的“灰箱”结构对于SISO 随机系统，被辨识模型()G z 为：12121212()()()1n n nn b z b z b z y z G z u z a z a z a z ------+++==++++ 其相应的差分方程为11()()()n ni i i i y k a y k i b u k i ===--+-∑∑若考虑被辨识系统或观测信息中含有噪声，被辨识模型可改写为11()()()()n ni i i i z k a y k i b u k i v k ===--+-+∑∑式中，z(k)为系统输出量的第k 次观测值；y(k)为系统输出量的第k 次真值，y(k-1)为系统输出量的第k-1次真值，以此类推；u(k)为系统的第k 个输入值，u(k-1)为系统的第k-1个输入值；v(k)为均值为0的不相关随机噪声。
(完整word版)2003版系统辨识最小二乘法大作业
西北工业大学系统辩识大作业题目：最小二乘法系统辨识一、 问题重述：用递推最小二乘法、加权最小二乘法、遗忘因子法、增广最小二乘法、广义最小二乘法、辅助变量法辨识如下模型的参数离散化有z^4 - 3.935 z^3 + 5.806 z^2 - 3.807 z + 0.9362---------------------------------------------- =z^4 - 3.808 z^3 + 5.434 z^2 - 3.445 z + 0.8187噪声的成形滤波器离散化有4.004e-010 z^3 + 4.232e-009 z^2 + 4.066e-009 z + 3.551e-010-----------------------------------------------------------------------------=z^4 - 3.808 z^3 + 5.434 z^2 - 3.445 z + 0.8187采样时间0.01s要求：1.用Matlab 写出程序代码；2.画出实际模型和辨识得到模型的误差曲线；3.画出递推算法迭代时各辨识参数的变化曲线；最小二乘法:在系统辨识领域中 ,最小二乘法是一种得到广泛应用的估计方法 ,可用于动态 ,静态 , 线性 ,非线性系统。
在使用最小二乘法进行参数估计时 ,为了实现实时控制 ,必须优化成参数递推算法 ,即最小二乘递推算法。
这种辨识方法主要用于在线辨识。
MATLAB 是一套高性能数字计算和可视化软件 ,它集成概念设计 ,算法开发 ,建模仿真 ,实时实现于一体 ,构成了一个使用方便、界面友好的用户环境 ,其强大的扩展功能为各领域的应用提供了基础。
对于一个简单的系统 ,可以通过分析其过程的运动规律 ,应用一些已知的定理和原理,建立数学模型 ,即所谓的“白箱建模 ”。
但对于比较复杂的生产过程 ,该建模方法有很大的局限性。
由于过程的输入输出信号一般总是可以测量的 ,而且过程的动态特性必然表现在这些输入输出数据中 ,那么就可以利用输入输出数据所提供的信息来建立过程的数学模型。
系统辨识matlab最小二乘法
一、 实验题目：最小二乘法在系统辨识中的应用二、 实验目的1.掌握系统辨识的理论、方法及应用2.熟练Matlab 下最小二乘法编程3.掌握M 序列产生方法三、 实验设备1、硬件设备：计算机配置，P4、32位CPU 、512M 内存 2、 软件设备： windows xp 操作系统 、matlab6.5软件包四、 实验原理最小二乘理论是有高斯（K.F.Gauss ）在1795年提出：“未知量的最大可能值是这样一个数值，它使各次实际观测值和计算值之间的差值的平方乘以度量其精度的数值以后的和最小。
”。
单输入单输出离散时间动态系统差分方程为：)()()()k (1i 1i k e i k u i k Z Z bni na i b a +-=-+∑∑==其中Z （k ）为输出变量，u(k)为输入变量，e(k)为偏差。
上式可以表示为)()()(-)k (1i 1i k e i k u i k Z Z bni na i b a +--=∑∑==各参数用矩阵表示 ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡--------------=)()1()()1()2()1()2()1()1()0()1()0(a a b a b a n l u l u n l z l z n u u n z z n u u n z z H （1） T l z z z Z )](),2(),1([，⋯= （2）其中l 为所需要采集的点数。
[]Tn e E )(........).........0(e = （3） []T bn anb a ........1.........1=θ （4）Z=H*θ+E ，E=Z-H*θ，根据最小二乘理论E 必须最小对上式进行求导，推出 Z H H H T T 1)(-=θ根据表达式Z H H H T T 1)(-=θ带入（1）（2）（4）即可求出a1....an b1.......bn 。
五、实验代码以及实验结果m=20; %置M序列总长度y1=1;y2=1;y3=1;y4=0;for i=1:mx1=xor(y3,y4);%异或运算x2=y1;x3=y2;x4=y3;if y4==0;u(i)=1;elseu(i)=-1;endy1=x1;y2=x2;y3=x3;y4=x4;endz=zeros(21,1);%定义输出观测值的长度21行*1列的0矩阵ZL=zeros(19,1);%定义输出观测值的长度19行*1列的0矩阵for k=3:21z(k)=-1.5*z(k-1)-z(k-2)+u(k-1)+3*u(k-2) ;%用理想输出值作为观测值ZL(k-2)=z(k);end%subplot(3,1,1) %画三行一列图形窗口中的第一个图形%stem(u) %画出输入信号u的图形%subplot(3,1,2) %画三行一列图形窗口中的第二个图形%i=1:1:16; %横坐标范围是1到16，步长为1%plot(i,z) %图形的横坐标是采样时刻i, 纵坐标是输出观测值z, 图形格式为连续曲线%subplot(3,1,3) %画三行一列图形窗口中的第三个图形%stem(z),grid %画出输出观测值z的图形，并显示坐标网格for m=2:20HL=[-z(m) -z(m-1) u(m) u(m-1)];for n=1:4;H((m-1),n)=HL(n);endendc1=H'*H;c2=inv(c1);c3=H'*ZL;c=c2*c3a1=c(1), a2=c(2), b1=c(3), b2=c(4)实验结果：c =1.50001.00001.00003.0000a1 =1.5000a2 =1.0000b1 =1b2 =3.0000六、实验结果分析通过实验结果可知所得的实验结果与待辨识的系统传递函数的系数很接近了。
MIMO系统最小二乘法参数辨识
转化为最小二乘形式
转换为标准最小二乘写法
式中，
转换最小二乘形式
2. 根据最小二乘法，求出对应的第i行对应的元素 3. 依次对每个输出yi，重复1-2 4. 将结果组合，得到相应的A和R 5. 由A和R，反求原来形式的矩阵，得到解。
最小二乘法只能处理如下的形式 其中，Y为Nx1的列向量，H为矩阵，r为参数列向量
最小二乘法特点
因为最小二乘发只能处理写成最小二乘表达式，
因此，最小二乘法很适合处理SISO的辨识。过程如下： 1. 将SISO系统写成差分方程表达形式，即连续传递函数离散化
最小二乘法
最小二乘法只能处理最小二乘形式，对系统辨识而言，只能直接处 理如下格式
s^2 + 2.545 s + 3.24 为了计算机便于计算，需要离散化后，得到离散格式的传递函数 dsys，对得到的dsys使用扫频信号u激励，得到输出y。根据得到 的u和y，辨识得到离散的传递函数dsysL，然后得到连续传递函数 sysL。
过程数据
辨识得到的结果
响应曲线验证
最小二乘法特点
最小二乘法特点
2. 上式中生成的a_i，b_i当做未知参数theta的分量 3. 写成最小二乘表达法
4. 求出a_i,b_i,得到离散传递函数（1步） 5. 由离散传递函数求得连续传递函数。
MIMO特点
状态空间表达式一般为
要想用最小二乘法辨识，有如下两个问题 1. 状态变量X需要消去
MIMO 最小二乘辨识
李曹磊
辨识目的
某些小型飞机不知道转动惯量、气动导数，又无法由传统方法获得 （或经费所限），可以通过辨识来求。从而为飞控系统设计提供自 然飞机模型。 对于某些飞机，可以验证传统方法得到的气动导数。
昆明理工数值分析大作业最小二乘法
昆明理工数值分析大作业最小二乘法最小二乘法（Least Squares Method）是数值分析中的一种重要方法，用于处理数据拟合问题。
在大作业中，我们将通过使用最小二乘法来拟合给定的数据，并解释其原理和应用。
最小二乘法是一种用于找到使得拟合曲线与数据点之间的误差最小化的方法。
使用最小二乘法进行数据拟合的基本思想是，找到一个函数，可以描述数据点的分布，并通过优化算法调整函数的参数，使得函数的拟合曲线与数据点的残差最小。
首先，我们需要确定拟合函数的形式。
在拟合直线的情况下，我们选择一条直线的方程 y = mx + b，其中 m 和 b 是需要衡量和优化的参数。
在更复杂的情况下，比如多项式拟合，拟合函数可以是二次函数、三次函数等。
最小二乘法的关键是定义误差函数或损失函数。
通常，最小二乘法使用残差平方和来作为误差函数。
残差是指拟合曲线与实际数据点之间的垂直距离。
对于一条直线来说，残差可以通过计算每个数据点在垂直方向上的距离来得到。
如果我们有n个数据点，那么残差平方和可以通过以下公式计算：S = Σ(yᵢ - (mxᵢ + b))²其中，(xᵢ,yᵢ)表示第i个数据点的坐标。
我们的目标是找到最佳的参数m和b，使得S最小化。
为了找到最小化残差平方和的解，可以使用最优化算法，如梯度下降法、牛顿法等。
这些算法根据误差函数的梯度（导数）来更新参数的值，直到达到最小化误差的目标。
最小二乘法在实际应用中有广泛的用途。
例如，在回归分析中，可以使用最小二乘法进行线性回归，以确定自变量和因变量之间的关系。
此外，最小二乘法还可以用于曲线拟合、信号处理、图像处理等领域。
在大作业中，你可以选择一个合适的数据集，并使用最小二乘法进行拟合。
你可以尝试不同的拟合函数和最优化算法，比较它们的性能和误差。
此外，你还可以进一步探索最小二乘法的应用领域，并说明其优缺点。
总之，最小二乘法是一种重要的数值分析方法，用于拟合数据并优化参数。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
系统辨识最小二乘法大作业系统辨识大作业最小二乘法及其相关估值方法应用学院:自动化学院专业:信息工程学号:2007302171姓名:马志强日期:2010.11.14基于最小二乘法的多种系统辨识方法研究 1. 最小二乘法的引出在系统辨识中用得最广泛的估计方法是最小二乘法(LS)。
设单输入-单输出线性定长系统的差分方程为(5.1.1)式中:为随机干扰;为理论上的输出值。
只有通过观测才能得到，在观测过程中往往附加有随机干扰。
的观测值可表示为(5.1.2)式中:为随机干扰。
由式(5.1.2)得(5.1.3)将式(5.1.3)带入式(5.1.1)得(5.1.4)我们可能不知道的统计特性，在这种情况下，往往把看做均值为0的白噪声。
设(5.1.5)则式(5.1.4)可写成(5.1.6)在观测时也有测量误差，系统内部也可能有噪声，应当考虑它们的影响。
因此假定不仅包含了的测量误差，而且包含了的测量误差和系统内部噪声。
假定是不相关随机序列(实际上是相关随机序列)。
现分别测出个随机输入值，则可写成个方程，即上述个方程可写成向量-矩阵形式(5.1.7)设则式(5.1.7)可写为(5.1.8)式中:为维输出向量;为维噪声向量;为维参数向量;为测量矩阵。
因此式(5.1.8)是一个含有个未知参数，由个方程组成的联立方程组。
如果，方程数少于未知数数目，则方程组的解是不定的，不能唯一地确定参数向量。
如果，方程组正好与未知数数目相等，当噪声时，就能准确地解出(5.1.9) 如果噪声，则(5.1.10) 从上式可以看出噪声对参数估计是有影响的，为了尽量较小噪声对估值的影响。
在给定输出向量和测量矩阵的条件下求系统参数的估值，这就是系统辨识问题。
可用最小二乘法来求的估值，以下讨论最小二乘法估计。
2. 最小二乘法估计算法设表示的最优估值，表示的最优估值，则有(5.1.11)写出式(5.1.11)的某一行，则有(5.1.12)设表示与之差，即-(5.1.13) 式中成为残差。
把分别代入式(5.1.13)可得残差。
设则有(5.1.14) 最小二乘估计要求残差的平方和为最小，即按照指数函数(5.1.15) 为最小来确定估值。
求对的偏导数并令其等于0可得(5.1.16)(5.1.17) 由式(5.1.17)可得的最小二乘估计(5.1.18) 3.递推最小二乘法为了实现实时控制，必须采用递推算法，这种辨识方法主要用于在线辨识。
设已获得的观测数据长度为，将式(5.1.8)中的和分别用来代替，即(5.3.1)用的最小二乘估计，则(5.3.2)设(5.3.5) 于是(5.3.6)如果再获得1组新的观测值和，则又增加1个方程(5.3.7)式中将式(5.3.1)和式(5.3.7)合并，并写成分块矩阵形式，可得(5.3.8) 根据上式可得到新的参数估值(5.3.9)式中根据矩阵求逆引理可以求得递推最小二乘法辨识公式(5.3.19)(5.3.20)(5.3.21)由于进行递推计算需要给出和初值和，通过计算证明，可以取初值:，，c是充分大的常数，为单位矩阵，则经过若干次递推之后能够得到较好的参数估计。
3. 辅助变量法辅助变量法是一种可克服最小二乘有偏估计的一种方法，对于原辨识方程(5.4.1)当是不相关随机序列时，最小二乘法可以得到参数向量的一致无偏估计。
但是，在实际应用中往往是相关随机序列。
假定存在着一个的矩阵满足约束条件(5.4.2)式中是非奇异的。
用乘以式(5.4.1)等号两边得(5.4.3) 由上式得(5.4.4) 如果取(5.4.5) 作为估值，则称为辅助变量估值，矩阵成为辅助变量矩阵，中的元素称为辅助变量。
常用的辅助变量法有递推辅助变量参数估计法，自适应滤波法，纯滞后等。
4. 广义最小二乘法广义最小二乘法是能克服最小二乘法有偏估计的另一种方法，这种方法计算比较复杂但效果比较好。
下面直接介绍广义最小二乘法的计算步骤:(1)应用得到的输入和输出数据和，按模型求出的最小二乘估计(2)计算残差(3)用残差代替，计算(4)计算和(5)应用得到的和按模型用最小二乘法重新估计，得到的第2次估值。
然后按步骤(2)计算残差，按步骤(3)重新估计，得到估值。
再按照步骤(4)计算和，按照步骤(5)求的第3次估值。
重复上述循环，之道的估值收敛为止。
5. 一种交替的广义最小二乘法求解技术(夏式法)这种方法是夏天长提出来的，又称夏式法。
以上讨论过的广义最小二乘法的特点在于系统的输入和输出信号反复过滤。
一下介绍的夏式法是一种交替的广义最小二乘法求解技术，它不需要数据反复过滤，因而计算效率较高。
这种方法可消去最小二乘估计中的偏差，而且由这种方法导出的计算方法也比较简单。
基于以上的几种方法，有(5.7.1) 因而有(5.7.2)应用最小二乘法可得到参数估值(5.7.3)可以推出(5.7.11)上式中的第1项是最小二乘估计，第2项是偏差项，所以必须准确计算。
为了准确计算，可采用迭代的方法。
6. 专题解答设但输入-单输出系统的差分方程为取真实值，输入数据如下所示k u(k) k u(k) k u(k)1 1.147 11 -0.958 21 0.4852 0.201 12 0.810 22 1.6333 -0.787 13 -0.044 23 0.0434 -1.159 14 0.947 24 1.3265 -1.052 15 -1.474 25 1.7066 0.866 16 -0.719 26 -0.3407 1.152 17 -0.086 27 0.8908 1.573 18 -1.099 28 1.1449 0.626 19 1.450 29 1.17710 0.433 20 1.151 30 -0.390用的真实值利用查分方程求出作为测量值，为均值为0，方差为0.1，0.5的不相关随机序列。
(1) 用最小二乘法估计参数。
(2) 用递推最小二乘法估计。
(3) 用辅助变量法估计参数。
(4) 设，用广义最小二乘法估计参数。
(5) 用夏式法估计参数。
(6) 详细分析和比较所获得的参数辨识结果，并说明上述参数便是方法的优缺点。
根据题目要求的解法，利用Matlab编程实现系统辨识的估值利用最小二乘法估计的结果如下:最小二乘法方差0.0001 1.6280 0.7028 0.3971 3.4491..6316 0.7059 0.3947 0.34941.6354 0.7120 0.3918 0.34631.6362 0.7082 0.3970 0.35271.6360 0.7165 0.3906 0.35341.6289 0.7046 0.3908 0.34370.001 1.1543 0.6766 0.4064 0.33041.5577 0.6371 0.3868 0.32491.6050 0.6860 0.3737 0.32441.6060 0.6816 0.3583 0.31671.6195 0.7030 0.3907 0.33661.5670 0.6572 0.3752 0.31400.01 1.3538 0.5010 0.3486 0.17090.8956 0.1637 0.4237 0.06971.0008 0.2036 0.4236 0.12681.3403 0.4707 0.3826 0.26151.0574 0.2289 0.3682 0.13171.1231 0.2963 0.3592 0.1506 0.1 1.1424 0.2710 0.3284 0.2216 1.0255 0.1736 0.3766 0.18440.8896 0.1109 0.3893 0.08130.8182 0.1114 0.4298 0.09230.8100 0.0153 0.4122 0.13520.7715 0.1311 0.4714 0.0640 0.5 0.9751 0.1017 0.0271 0.0792 0.8938 0.0740 0.3484 0.16340.1927 0.0197 0.3762 0.15210.5506 0.0392 0.6510 0.05840.7560 0.0494 0.3372 0.14370.9459 0.1377 0.3815 0.1853部分程序运行结果递推最小二乘法方差 0.0001 1.6754 0.6787 0.5207 0.39241.3076 0.2900 0.0075 0.32841.5146 0.6963 1.1401 0.16391.5733 0.7782 0.4149 0.71101.1602 0.4753 0.6736 0.31591.2091 0.3192 0.5277 0.0275 0.001 1.4767 0.4040 0.2679 0.5512 1.6259 0.7594 0.3253 0.37801.5393 0.4757 0.1268 0.43461.1548 0.1700 0.1926 0.82510.8858 0.0760 0.3385 0.04061.4129 0.3127 0.0992 0.83800.01 1.3485 0.3445 0.3194 0.37101.1639 0.3296 0.7813 0.21621.9946 1.2323 1.4852 0.03041.3924 0.3543 0.3319 0.45721.3982 0.3608 0.7773 0.31521.6346 0.7229 0.5780 0.3947 0.1 1.5624 0.7132 0.4422 0.4112 1.7335 0.7152 0.0844 0.63991.4763 0.5366 0.3255 0.31161.4477 0.3489 0.2218 0.22651.6216 0.7082 0.6595 0.42751.5105 0.4000 0.0113 0.2213 0.5 1.7927 0.9411 0.2730 0.3471 1.5556 0.8877 0.5972 0.12171.7868 1.2538 1.1248 0.21001.5733 0.7434 0.3589 0.13871.3193 0.6084 1.2971 0.30291.5959 0.5386 0.0141 0.6947部分程序运行结果:辅助变量法方差 0.0001 1.7799 0.8588 0.4147 0.40461.3076 0.2900 0.0075 0.32841.6735 0.7578 0.4060 0.34841.5812 0.6546 0.3771 0.36111.6657 0.7469 0.3772 0.35611.5281 0.6509 0.3645 0.3302 0.001 1.6295 0.6775 0.4082 0.3501 1.6425 0.7305 0.3937 0.35431.5595 0.6052 0.3563 0.33621.4145 0.4925 0.4021 0.28731.6371 0.7270 0.3681 0.34181.3539 0.4733 0.3906 0.24890.01 1.3451 0.4700 0.3822 0.26801.3657 0.4893 0.4349 0.24961.3702 0.5009 0.4388 0.26111.1884 0.3707 0.3457 0.15211.3636 0.5330 0.4135 0.20641.3158 0.4836 0.4562 0.2169 0.1 1.5545 0.6167 0.4104 0.3512 1.5900 0.6648 0.4098 0.36641.6610 0.7029 0.4001 0.34951.5104 0.6014 0.3968 0.32561.5620 0.6496 0.3950 0.30971.4418 0.5706 0.4183 0.2783 0.5 1.4952 0.5704 0.3769 0.3483 1.5592 0.6541 0.4387 0.33301.3637 0.5302 0.3865 0.19161.5543 0.6324 0.3208 0.29291.5385 0.5994 0.3661 0.35761.3511 0.4839 0.3412 0.2208部分程序运行结果:广义最小二乘法方差 0.0001 1.6451 0.7182 0.3895 0.35031.6468 0.7205 0.3866 0.34921.6308 0.7038 0.3852 0.34571.6467 0.7188 0.3898 0.35441.6411 0.7156 0.3898 0.34771.6441 0.7185 0.3896 0.3488 0.001 1.7027 0.7678 0.4059 0.3761 1.6477 0.7166 0.3973 0.36431.6391 0.7113 0.3836 0.35041.6586 0.7166 0.3791 0.36721.6831 0.7472 0.3994 0.37311.5964 0.6850 0.3655 0.32310.01 1.6826 0.7469 0.3859 0.36001.7245 0.7696 0.3565 0.36421.6577 0.7186 0.3910 0.38201.6656 0.7263 0.3425 0.34131.6903 0.7392 0.3633 0.38381.6942 0.7439 0.3834 0.3880 0.1 0.8996 0.1515 0.4435 0.0604 1.0380 0.3044 0.3410 0.04641.6571 0.6569 0.4589 0.42591.2180 0.4584 0.3152 0.10300.1957 0.0974 0.2164 0.10351.6859 0.7232 0.4144 0.4404 0.5 1.4053 0.6703 0.3672 0.1372 1.5783 0.6148 0.3512 0.41731.0233 0.0658 0.2035 0.40601.5973 0.8103 0.3544 0.21311.6191 0.6709 0.3111 0.48661.0804 0.2210 0.7500 0.5185部分程序运行结果:夏式法方差1.6367 0.7090 0.3910 0.3511 0.00011.6451 0.7184 0.3876 0.34721.3618 0.2159 0.8205 0.58511.6291 0.7075 0.3911 0.34641.6173 0.6929 0.3955 0.34391.6307 0.7066 0.3912 0.3449 0.001 1.5314 0.6407 0.4016 0.3092 1.5632 0.6575 0.4204 0.33191.6172 0.6946 0.3925 0.32581.5670 0.6489 0.4016 0.32331.5326 0.6243 0.3995 0.30701.6251 0.6946 0.3621 0.33600.01 0.8815 0.1089 0.3734 0.04401.3763 0.4972 0.4121 0.26541.5027 0.5952 0.3854 0.32771.3421 0.4268 0.3978 0.27901.2454 0.3024 0.3669 0.25251.3361 0.4741 0.4999 0.29000.1 0.6639 0.1074 0.6847 0.08271.1406 0.3154 0.4259 0.23670.1791 0.2916 0.9454 0.21261.2068 0.3827 0.7460 0.32200.6412 0.0704 0.5802 0.08771.0384 0.2878 0.7493 0.37930.9330 0.0752 0.9741 0.5882 0.50.6291 0.1686 0.4679 0.78690.5986 0.2792 0.2605 0.20080.4441 0.8033 0.7442 1.76670.6216 0.1440 0.5821 0.49060.6576 0.0999 1.0002 0.0263部分程序运行结果:结论:通过编程计算，获得在噪声方差比较小的情况下，各种方法所获得的估值比较理想，但随着噪声方差的增大，估值的偏差随之增大，横向比较看来夏式法与广义最小二乘法能够更好地还原参数值，当观测值足够多时，各种方法都能很好地反映参数真实值。