最小二乘法MATLAB程序及结果
最小二乘法曲线拟合的Matlab程序

方便大家使用的最小二乘法曲线拟合的Matlab程序非常方便用户使用,直接按提示操作即可;这里我演示一个例子:(红色部分为用户输入部分,其余为程序运行的结果,结果图为Untitled.fig,Untitled2.fig) 请以向量的形式输入x,y.x=[1,2,3,4]y=[3,4,5,6]通过下面的交互式图形,你可以事先估计一下你要拟合的多项式的阶数,方便下面的计算.polytool()是交互式函数,在图形上方[Degree]框中输入阶数,右击左下角的[Export]输出图形回车打开polytool交互式界面回车继续进行拟合输入多项式拟合的阶数m = 4Warning: Polynomial is not unique; degree >= number of data points. > In polyfit at 72In zxecf at 64输出多项式的各项系数a = 0.0200000000000001a = -0.2000000000000008a = 0.7000000000000022a = 0.0000000000000000a = 2.4799999999999973输出多项式的有关信息 SR: [4x5 double]df: 0normr: 2.3915e-015Warning: Zero degrees of freedom implies infinite error bounds.> In polyval at 104In polyconf at 92In zxecf at 69观测数据拟合数据x y yh1.0000 3.0000 3.00002.0000 4.0000 4.00003 5 54.0000 6.0000 6.0000剩余平方和 Q = 0.000000标准误差 Sigma = 0.000000相关指数 RR = 1.000000请输入你所需要拟合的数据点,若没有请按回车键结束程序.输入插值点x0 = 3输出插值点拟合函数值 y0 = 5.0000>>结果:untitled.figuntitled2.fig一些matlab优化算法代码的分享代码的目录如下:欢迎讨论1.约束优化问题:minRosen(Rosen梯度法求解约束多维函数的极值)(算法还有bug) minPF(外点罚函数法解线性等式约束)minGeneralPF(外点罚函数法解一般等式约束)minNF(内点罚函数法)minMixFun(混合罚函数法)minJSMixFun(混合罚函数加速法)minFactor(乘子法)minconPS(坐标轮换法)(算法还有bug)minconSimpSearch(复合形法)2.非线性最小二乘优化问题minMGN(修正G-N法)3.线性规划:CmpSimpleMthd(完整单纯形法)4.整数规划(含0-1规划)DividePlane(割平面法)ZeroOneprog(枚举法)5.二次规划QuadLagR(拉格朗日法)ActivedeSet(起作用集法)6.辅助函数(在一些函数中会调用)minNT(牛顿法求多元函数的极值)Funval(求目标函数的值)minMNT(修正的牛顿法求多元函数极值)minHJ(黄金分割法求一维函数的极值)7.高级优化算法1)粒子群优化算法(求解无约束优化问题)1>PSO(基本粒子群算法)2>YSPSO(待压缩因子的粒子群算法)3>LinWPSO(线性递减权重粒子群优化算法)4>SAPSO(自适应权重粒子群优化算法)5>RandWSPO(随机权重粒子群优化算法)6>LnCPSO(同步变化的学习因子)7>AsyLnCPSO(异步变化的学习因子)(算法还有bug)8>SecPSO(用二阶粒子群优化算法求解无约束优化问题)9>SecVibratPSO(用二阶振荡粒子群优化算法求解五约束优化问题)10>CLSPSO(用混沌群粒子优化算法求解无约束优化问题)11>SelPSO(基于选择的粒子群优化算法)12>BreedPSO(基于交叉遗传的粒子群优化算法)13>SimuAPSO(基于模拟退火的粒子群优化算法)2)遗传算法1>myGA(基本遗传算法解决一维约束规划问题)2>SBOGA(顺序选择遗传算法求解一维无约束优化问题)3>NormFitGA(动态线性标定适应值的遗传算法求解一维无约束优化问题)4>GMGA(大变异遗传算法求解一维无约束优化问题)5>AdapGA(自适应遗传算法求解一维无约束优化问题)6>DblGEGA(双切点遗传算法求解一维无约束优化问题)7>MMAdapGA(多变异位自适应遗传算法求解一维无约束优化问题)自己编写的马尔科夫链程序A 代表一组数据序列一维数组本程序的操作对象也是如此t=length(A); % 计算序列“A”的总状态数B=unique(A); % 序列“A”的独立状态数顺序,“E”E=sort(B,'ascend');a=0;b=0;c=0;d=0;for j=1:1:ttLocalization=find(A==E(j)); % 序列“A”中找到其独立状态“E”的位置for i=1:1:length(Localization)if Localization(i)+1>tbreak; % 范围限定elseif A(Localization(i)+1)== E(1)a=a+1;elseif A(Localization(i)+1)== E(2)b=b+1;elseif A(Localization(i)+1)== E(3)c=c+1;% 依此类推,取决于独立状态“E”的个数elsed=d+1;endendT(j,1:tt)=[a,b,c,d]; % “T”为占位矩阵endTT=T;for u=2:1:ttTT(u,:)= T(u,:)- T(u-1,:);endTT; % 至此,得到转移频数矩阵Y=sum(TT,2);for uu=1:1:ttTR(uu,:)= TT(uu,:)./Y(uu,1);endTR % 最终得到马尔科夫转移频率/概率矩阵% 观测序列马尔科夫性质的检验:N=numel(TT);uuu=1;Col=sum(TT,2); % 对列求和Row=sum(TT,1); % 对行求和Total=sum(Row); % 频数总和for i=1:1:ttfor j=1:1:ttxx(uuu,1)=sum((TT(i,j)-(Row(i)*Col(j))./Total).^2./( (Row(i)*Col(j)). /Total));uuu=uuu+1; % 计算统计量x2endendxx=sum(xx)。
matlab多元一次方程最小二乘拟合求取系数

matlab多元一次方程最小二乘拟合求取系数最小二乘法是一种数学方法,可用于数据拟合以及误差分析。
在科学与工程中,我们经常需要解决数据拟合的问题,而最小二乘法便是其中最常用的方法之一。
本文将介绍如何使用MATLAB进行多元一次方程最小二乘拟合并求取系数。
第一步:准备数据首先,我们需要准备一组数据以进行最小二乘拟合。
这组数据需要是多个变量之间的关系,并且这些变量需要满足线性关系。
对于这个问题,我们可以先简化成一个二元一次方程y=ax+b。
这个方程可以表示成矩阵形式:```Y = [y1;y2;...;ym]X = [1,x1;1,x2;...;1,xm]B = [b;a]```其中,Y是一个m行一列的向量,表示对应的y值;X是一个m 行两列的矩阵,第一列为1表示截距项,第二列为x值;B是一个两行一列的向量,表示最终的系数。
第二步:计算最小二乘法接下来,我们需要使用MATLAB求解这个最小二乘问题。
我们可以使用MATLAB内置的regress函数,它可以帮助我们求解系数B。
具体使用方法如下:```B = regress(Y,X)```这个函数会返回一个2行1列的向量,也就是系数B的值。
第三步:验证结果最后,我们需要验证我们拟合出的结果是否可靠。
我们可以使用拟合残差来评估我们的拟合效果,同时也可以用图形的方式来直观地观察。
对于残差,可以使用如下代码来计算:```e = Y - X * B```这个代码会计算出每一行数据的残差。
我们可以使用hist函数来显示残差分布的直方图:```hist(e)```对于图形,我们可以使用plot函数来绘制拟合曲线。
具体代码如下:```plot(X(:,2),Y,'o')hold onplot(X(:,2),X * B,'-')hold off```这个代码会在同一个坐标系内绘制出数据点以及拟合曲线,以直观地观察拟合效果。
以上就是在MATLAB中进行多元一次方程最小二乘拟合的步骤。
最小二乘法matlab程序

最小二乘法(Least Squares Method,LSM)是一种数值计算方法,用于拟合曲线,求解未知参数的值。
它的基本思想是,通过求解最小二乘误差的最优解,来拟合曲线,从而求得未知参数的值。
本文将介绍最小二乘法在Matlab中的实现原理及程序编写。
一、最小二乘法的原理最小二乘法是一种数值计算方法,它的基本思想是,通过求解最小二乘误差的最优解,来拟合曲线,从而求得未知参数的值。
最小二乘法的基本原理是:给定一组数据点,用直线拟合这组数据点,使得拟合直线与这组数据点的误差的平方和最小。
具体地说,假设有一组数据点,其中每个数据点都可表示为(x_i, y_i),i=1,2,3,...,n,其中x_i和y_i分别表示第i个数据点的横纵坐标。
拟合这组数据点的直线通常用一元线性函数表示,即y=ax+b,其中a和b是未知参数。
最小二乘法的思想是:求出使误差的平方和最小的a和b,即求出最优解。
二、Matlab程序编写1. 准备工作首先,我们需要准备一组数据点,每个数据点都可表示为(x_i, y_i),i=1,2,3,...,n,其中x_i和y_i分别表示第i个数据点的横纵坐标。
例如,我们可以准备一组数据点:x=[1,2,3,4,5];y=[2,4,6,8,10];2. 程序编写接下来,我们就可以开始编写Matlab程序了。
首先,我们需要定义一个一元线性函数,用于拟合这组数据点。
函数的形式为:y=ax+b,其中a和b是未知参数。
%定义函数f=@(a,b,x)a*x+b;然后,我们需要定义一个误差函数,用于计算拟合直线与这组数据点的误差的平方和。
%定义误差函数error=@(a,b)sum((y-f(a,b,x)).^2);最后,我们就可以使用Matlab提供的fminsearch函数,求解最小二乘误差的最优解,即求出最优a和b的值。
%求解最优解[a,b]=fminsearch(error,[1,1]);经过上面的程序编写,我们就可以求得未知参数a和b的最优值。
各种最小二乘法汇总(算例及MATLAB程序)

u(400)
⎟ ⎠
Matlab程序见附录 1。
1.2. 递推最小二乘算法
递推最小二乘算法公式:
^
^
^
θ (k) = θ (k −1) + K (k)[z(k) − h' (k)θ (k −1)]
K (k) = P(k −1)h(k)[h' (k)P(k −1)h(k) + 1 ]−1
(1.2)
Λ(k)
b2 a1 50 100 150 200 250 300 350 400 450
图 1 一般最小二乘参数过渡过程
4
盛晓婷 最小二乘算法总结报告
估计方差变化过程
100
90
80
70
60
50
40
30
20
10
0
0
50 100 150 200 250 300 350 400 450
图 2 一般最小二乘方差变化过程
1.1. 一次计算最小二乘算法
⎛ ⎜
^
a1
⎞ ⎟
⎛ -1.4916 ⎞
^
θ
LS
=
⎜^ ⎜ a2 ⎜^ ⎜ b1
⎟ ⎟ ⎟ ⎟
=
(
H
T L
H
L
)−1
H
T L
Z
L
=
⎜ ⎜
0.7005
⎟ ⎟
⎜1.0364 ⎟
⎜
⎟
(1.1)
⎜⎜⎝
^
b2
⎟⎟⎠
⎝ 0.4268 ⎠
⎛ Z (3) ⎞
⎛ hT (3) ⎞ ⎛ −Z (2) −Z (1)
图 1 一般最小二乘参数过渡过程 .....................................................4 图 2 一般最小二乘方差变化过程 ....................................................5 图 3 遗忘因子法参数过渡过程 ........................................................7 图 4 遗忘因子法方差变化过程 ........................................................8 图 5 限定记忆法参数过渡过程 ......................................................10 图 6 限定记忆法方差变化过程 ......................................................10 图 7 偏差补偿最小二乘参数过渡过程 ..........................................12 图 8 偏差补偿最小二乘方差变化过程 ..........................................12 图 9 增广最小二乘辨识模型 ..........................................................13 图 10 增广最小二乘参数过渡过程 ................................................14 图 11 广义最小二乘参数过渡过程 ................................................16 图 12 广义最小二乘方差变化过程 ................................................16 图 13 辅助变量法参数过渡过程 ....................................................18 图 14 辅助变量法方差变化过程 ....................................................18 图 15 二步法参数过渡过程 ............................................................20 图 16 二步法方差变化过程 ............................................................20
matlab最小二乘解方程

matlab最小二乘解方程最小二乘法是求解线性方程组的一种有效方法,可以通过最小化误差平方和来得到最优解。
在MATLAB中,我们可以使用“\”操作符或者使用“pinv”函数来求解一个线性方程组的最小二乘解。
以下是关于如何在MATLAB中使用最小二乘法来求解线性方程组的详细内容:1. 使用“\”操作符使用“\”操作符可以很方便地求解一个线性方程组的最小二乘解。
例如,假设我们有一个由n个方程组成的线性方程组:Ax = b其中,A是一个m ×n的矩阵,x是一个n维向量,b是一个m维向量。
则它的最小二乘解为:x = (A' A)^(-1) A' b在MATLAB中,我们可以通过以下代码实现最小二乘解:A = [1 1 1; 2 3 4; 4 5 7; 5 6 8];b = [1; 2; 3; 4];x = A \ b;其中,反斜杠符号“\”表示求解线性方程组的最小二乘解。
2. 使用“pinv”函数除了使用“\”操作符,我们也可以使用MATLAB中的“pinv”函数来求解一个线性方程组的最小二乘解。
例如,我们可以通过以下代码实现最小二乘解:A = [1 1 1; 2 3 4; 4 5 7; 5 6 8];b = [1; 2; 3; 4];x = pinv(A) * b;其中,pinv函数表示求矩阵A的伪逆矩阵。
使用“pinv”函数来求解线性方程组的最小二乘解与使用“\”操作符的结果是等价的。
需要注意的是,在使用最小二乘法来求解线性方程组时,矩阵A的列应该是线性无关的,否则可能会出现唯一最小二乘解不存在的情况。
综上所述,MATLAB中使用最小二乘法来求解线性方程组非常简单。
我们可以通过“\”操作符或者“pinv”函数来求解一个线性方程组的最小二乘解。
最小二乘法matlab程序

最小二乘法matlab程序最小二乘法是一种统计模型,它可以被用来拟合一元函数数据,或者拟合非线性曲线。
它的基本思想是找到一组参数,使得拟合的曲线与实际数据的差距最小。
本文将介绍如何使用Matlab实现一个最小二乘法的程序,并与现有的一些现成的最小二乘法的matlab程序进行比较,找出其优缺点。
首先,要使用最小二乘法拟合曲线,需要准备一组输入数据,一般可以将其表示为两个向量,分别是自变量x和因变量y。
这些数据可以是由测量和实验得到的,也可以是由人工输入的,但无论如何都要确保它们的准确性。
接下来,就可以使用Matlab输入数据进行处理,用最小二乘法计算出最拟合的曲线及其参数。
具体步骤主要分为三步:第一步是计算输入数据的均值和方差,包括自变量x和因变量y的均值和方差;第二步是计算自变量x和因变量y的关系,即最小二乘拟合曲线的系数;第三步是验证拟合的曲线的准确性,如果不满意,可以重新调整参数,以获得较好的拟合效果。
此外,Matlab除了提供自带的最小二乘法函数外,还支持第三方开发者开发现成的matlab程序,用于解决最小二乘法的问题。
这些程序中有一些是开源的,另一些则是出售的。
其中开源的有LEAST,CURVEFIT,CURVEFITTOOL等,而出售的有MATLAB Curve Fitting Toolbox,Optimization Toolbox和Statistics Toolbox等。
它们的突出特点是速度快,代码简洁,容易上手,适用于多种拟合类型。
然而,各种matlab程序也有自身的缺点,最明显的就是当输入数据非常庞大时,它们的计算能力就无法跟上,速度就会变慢。
此外,使用出售的matlab程序可能相对昂贵,而且有时需要安装某些复杂的库文件,这也是一种麻烦。
因此,使用最小二乘法拟合曲线时,可以参考现有的matlab程序,也可以自己编写matlab代码,同时要考虑到程序的可靠性、效率和可行性。
本文介绍的matlab程序的最大优势是它不需要依赖第三方的软件,而且能够满足大多数用户的需求,使得最小二乘法可以在短时间内被成功运用。
matlab 最小二乘拟合直线并输出直线方程

在Matlab中,最小二乘法是一种常见的数学拟合技术,可以用来拟合直线,曲线甚至更复杂的函数。
通过最小二乘法,可以找到最适合数据点的直线方程,从而能够更好地分析和预测数据之间的关系。
在本文中,我将详细介绍如何在Matlab中使用最小二乘法来拟合直线,并输出直线方程。
我们需要准备一组数据点。
假设我们有一组横坐标和纵坐标的数据点,分别用变量x和y表示。
接下来,我们可以使用Matlab中的polyfit函数来进行最小二乘拟合。
该函数的语法如下:```matlabp = polyfit(x, y, 1);```其中,x和y分别代表数据点的横坐标和纵坐标,而1代表要拟合的直线的次数,即一次函数。
执行该语句后,变量p将会存储拟合出的直线的系数,即直线方程y = ax + b中的a和b。
在接下来的内容中,我将详细讨论如何通过最小二乘法拟合直线,并输出直线方程。
具体而言,我们将从如何准备数据、使用polyfit函数进行拟合、得到直线方程以及如何应用和解释直线拟合结果等方面进行全面分析。
一、数据准备在使用最小二乘法拟合直线之前,首先要准备一组数据点。
这些数据点应该是具有一定规律性的,从而能够通过直线拟合来揭示数据之间的关系。
在这一部分,我将详细介绍如何准备数据,并重点关注数据的合理性和可靠性。
1.1 数据收集要拟合直线,首先需要收集一组数据点。
这些数据点可以来源于实验观测、实际测量或者模拟计算等方式。
在收集数据时,需要保证数据的准确性和完整性。
还需要考虑数据的分布范围和密度,以便更好地反映数据之间的关系。
1.2 数据预处理在拟合直线之前,通常需要对数据进行一定的预处理。
这可能包括去除异常值、处理缺失数据,甚至进行数据变换等操作。
在这一步中,我将介绍如何进行数据预处理,并强调预处理对最终拟合结果的影响。
二、最小二乘拟合当数据准备工作完成后,就可以使用polyfit函数进行最小二乘拟合了。
在这一部分,我将详细介绍polyfit函数的使用方法,并解释其背后的数学原理。
matlab拟合最小二乘法

Ja31 = 49967500*a1 + 1089000*a2 + 25300*a3 + 660*a4 - 27131/100000
Ja41 = 1089000*a1 + 25300*a2 + 660*a3 + 24*a4 - 3987/500000 解线性方程组 Ja11 =0,Ja21 =0,Ja31 =0,Ja41 =0,输入下列程序 >> A=[117186437500, 2386725000,49967500,1089000; 2386725000, 49967500, 1089000,25300; 49967500,1089000,25300, 660; 1089000, 25300, 660,24]; B=[274377591928296252150123/576460752303423488000,31331074233255294718193/28823 03761517117440000,7819978335372091569501/28823037615171174400000, 2298349019433749545307/288230376151711744000000]; C=B/A, f=poly2sym(C) 运行即可得 C= 1.0e-004 * 0.0000 -0.0052 0.2634 0.0178
J= 58593218750*a1^2 + 2386725000*a1*a2 + 49967500*a1*a3 + 1089000*a1*a4 (274377591928296252150123*a1)/576460752303423488000 + 24983750*a2^2 + 1089000*a2*a3 + 25300*a2*a4 - (31331074233255294718193*a2)/2882303761517117440000 + 12650*a3^2 + 660*a3*a4 - (7819978335372091569501*a3)/28823037615171174400000 + 12*a4^2 (2298349019433749545307*a4)/288230376151711744000000 + 520374483464852566590953249225508026224249/33230699894622896822595176507008614 4000000000000 四、求 a1、a2、a3、a4 使 J 达到最小,分别对 a1、a2、a3、a4 求偏导数,使之等于 0 程序如下 syms a1 a2 a3 a4 J=58593218750*a1^2+2386725000*a1*a2+49967500*a1*a3+1089000*a1*a4-(2743775919282 96252150123*a1)/576460752303423488000+24983750*a2^2+1089000*a2*a3 + 25300*a2*a4 (31331074233255294718193*a2)/2882303761517117440000+12650*a3^2+660*a3*a4-(781997 8335372091569501*a3)/28823037615171174400000+12*a4^2-(2298349019433749545307*a4) /288230376151711744000000+520374483464852566590953249225508026224249/332306998 946228968225951765070086144000000000000 Ja1=diff(J,a1); Ja2=diff(J,a2); Ja3=diff(J,a3); Ja4=diff(J,a4); Ja11=simple(Ja1), Ja21=simple(Ja2), Ja31=simple(Ja3), Ja41=simple(Ja4) 运行得
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
最小二乘递推算法的MATLAB仿真针对辨识模型,有z(k)-+a1*z(k-1)+a2*z(k-2)=b1*u(k-1)+b2*u(k-2)+v(k)模型结构,对其进行最小二乘递推算法的MATLAB仿真,对比真值与估计值。
更改a1、a2、b1、b2参数,观察结果。
仿真对象:z(k)-1.5*z(k-1)+0.7*z(k-2)=u(k-1)+0.5*u(k-2)+v(k)程序如下:L=15;y1=1;y2=1;y3=1;y4=0; %四个移位寄存器的初始值for i=1:L; %移位循环x1=xor(y3,y4);x2=y1;x3=y2;x4=y3;y(i)=y4; %取出作为输出信号,即M序列if y(i)>0.5,u(i)=-0.03; %输入信号else u(i)=0.03;endy1=x1;y2=x2;y3=x3;y4=x4;endfigure(1);stem(u),grid onz(2)=0;z(1)=0;for k=3:15;z(k)=1.5*z(k-1)-0.7*z(k-2)+u(k-1)+0.5*u(k-2); %输出采样信号endc0=[0.001 0.001 0.001 0.001]'; %直接给出被识别参数的初始值p0=10^6*eye(4,4); %直接给出初始状态P0E=0.000000005;c=[c0,zeros(4,14)];e=zeros(4,15);for k=3:15; %开始求kh1=[-z(k-1),-z(k-2),u(k-1),u(k-2)]';x=h1'*p0*h1+1;x1=inv(x);k1=p0*h1*x1; %开始求k的值d1=z(k)-h1'*c0;c1=c0+k1*d1;e1=c1-c0;e2=e1./c0; %求参数的相对变化e(:,k)=e2;c0=c1;c(:,k)=c1;p1=p0-k1*k1'*[h1'*p0*h1+1]; %求出P(k)的值p0=p1;if e2<=E break;endendc,e %显示被辨识参数及其误差情况a1=c(1,:);a2=c(2,:);b1=c(3,:);b2=c(4,:);ea1=e(1,:);ea2=e(2,:);eb1=e(3,:);eb2=e(4,:);figure(2);i=1:15;plot(i,a1,'r',i,a2,':',i,b1,'g',i,b2,':')title('Parameter Identification with Recursive Least Squares Method')figure(3);i=1:15;plot(i,ea1,'r',i,ea2,'g',i,eb1,'b',i,eb2,'r:')title('Identification Precision')程序运行结果:p0 =1000000 0 0 00 1000000 0 00 0 1000000 00 0 0 1000000c =Columns 1 through 90.0010 0 0.0010 -0.4984 -1.2325 -1.4951 -1.4962 -1.4991 -1.49980.0001 0 0.0001 0.0001 -0.2358 0.6912 0.6941 0.6990 0.69980.0010 0 0.2509 1.2497 1.0665 1.0017 1.0020 1.0002 0.99990.0010 0 -0.2489 0.7500 0.5668 0.5020 0.5016 0.5008 0.5002Columns 10 through 15-1.4999 -1.5000 -1.5000 -1.5000 -1.4999 -1.49990.6999 0.7000 0.7000 0.7000 0.7000 0.70000.9998 0.9999 0.9999 0.9999 0.9999 0.99990.5002 0.5000 0.5000 0.5000 0.5000 0.5000e =1.0e+003 *Columns 1 through 90 0 0 -0.4994 0.0015 0.0002 0.0000 0.0000 0.00000 0 0 0 -2.3592 -0.0039 0.0000 0.0000 0.00000 0 0.2499 0.0040 -0.0001 -0.0001 0.0000 -0.0000 -0.00000 0 -0.2499 -0.0040 -0.0002 -0.0001 -0.0000 -0.0000 -0.0000Columns 10 through 150.0000 0.0000 0.0000 -0.0000 -0.0000 0.00000.0000 0.0000 -0.0000 0.0000 0.0000 0.0000-0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000-0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000程序运行曲线:图1.输入信号图2.a1,a2,b1,b2辨识仿真结果图3. a1,a2,b1,b2各次辨识结果收敛情况分析:由运行结果可看出,输出观测值没有任何噪声成分时,辨识结果最大相对误差达到3位数。
更改仿真对象参数:z(k)=z(k-1)-1.5*z(k-2)+u(k-1)+1.5*u(k-2);程序如下:L=15;y1=1;y2=1;y3=1;y4=0;for i=1:L;x1=xor(y3,y4);x2=y1;x3=y2;x4=y3;y(i)=y4;if y(i)>0.5,u(i)=-0.03;else u(i)=0.03;endy1=x1;y2=x2;y3=x3;y4=x4;endfigure(1);stem(u),grid onz(2)=0;z(1)=0;for k=3:15;z(k)=z(k-1)-1.5*z(k-2)+u(k-1)+1.5*u(k-2);endc0=[0.001 .0001 0.001 0.001]';p0=10^6*eye(4,4)E=0.000000005;c=[c0,zeros(4,14)];e=zeros(4,15);for k=3:15;h1=[-z(k-1),-z(k-2),u(k-1),u(k-2)]';x=h1'*p0*h1+1;x1=inv(x);k1=p0*h1*x1;d1=z(k)-h1'*c0;c1=c0+k1*d1;e1=c1-c0;e2=e1./c0;e(:,k)=e2;c0=c1;c(:,k)=c1;p1=p0-k1*k1'*[h1'*p0*h1+1];p0=p1;if e2<=E break;endendc,ea1=c(1,:);a2=c(2,:);b1=c(3,:);b2=c(4,:);ea1=e(1,:);ea2=e(2,:);eb1=e(3,:);eb2=e(4,:);figure(2);i=1:15;plot(i,a1,'r',i,a2,':',i,b1,'g',i,b2,':')title('Parameter Identification with Recursive Least Squares Method')figure(3);i=1:15;plot(i,ea1,'r',i,ea2,'g',i,eb1,'b',i,eb2,'r:')title('Identification Precision')程序运行结果:p0 =1000000 0 0 00 1000000 0 00 0 1000000 00 0 0 1000000c =Columns 1 through 90.0010 0 0.0010 0.4447 -1.2248 -0.9996 -0.9999 -1.0001 -1.00010.0001 0 0.0001 0.0001 0.3757 1.4987 1.4997 1.50001.50000.0010 0 -0.2489 0.6385 1.0560 1.0000 1.0001 0.99990.99990.0010 0 0.2509 1.1382 1.5557 1.4997 1.4995 1.49951.4996Columns 10 through 15-1.0001 -1.0000 -1.0000 -1.0000 -1.0000 -1.00001.5000 1.5000 1.5000 1.5000 1.5000 1.50000.9999 0.9999 0.9999 0.9999 0.9999 0.99991.4996 1.4996 1.4999 1.4999 1.4999 1.4999e =1.0e+003 *Columns 1 through 90 0 0 0.4437 -0.0038 -0.0002 0.0000 0.0000 -0.00000 0 0 0 3.7564 0.0030 0.0000 0.0000 -0.00000 0 -0.2499 -0.0036 0.0007 -0.0001 0.0000 -0.0000 0.00000 0 0.2499 0.0035 0.0004 -0.0000 -0.0000 0.0000 0.0000Columns 10 through 150.0000 -0.0000 -0.0000 0.0000 -0.0000 0.00000.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000-0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000-0.0000 -0.0000 0.0000 0.0000 0.0000 -0.0000程序运行曲线:图4. 输入信号图5. a1,a2,b1,b2辨识仿真结果图6. a1,a2,b1,b2各次辨识结果的收敛情况结果总结及问题分析:辨识结果与程序运行曲线表明,大约递推到五步时,参数辨识的结果基本达到稳定状态。
在参数更改前后运行结果可以看出,输出观测值在没有任何噪声成分下,估计值与真值最大相对误差达到3位数。