0ce3455c26a6147b657124b41cfe312dc0b44b7e
[scilab.git] / scilab / modules / signal_processing / help / en_US / identification / rpem.xml
1 <?xml version="1.0" encoding="UTF-8"?>
2 <refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:lang="en" xml:id="rpem">
3     <refnamediv>
4         <refname>rpem</refname>
5         <refpurpose>Recursive Prediction-Error Minimization estimation</refpurpose>
6     </refnamediv>
7     <refsynopsisdiv>
8         <title>Calling Sequence</title>
9         <synopsis>[w1,[v]]=rpem(w0,u0,y0,[lambda,[k,[c]]])</synopsis>
10     </refsynopsisdiv>
11     <refsection>
12         <title>Arguments</title>
13         <variablelist>
14             <varlistentry>
15                 <term>w0</term>
16                 <listitem>
17                     <para>
18                         <literal>list(theta,p,l,phi,psi)</literal> where:
19                     </para>
20                     <variablelist>
21                         <varlistentry>
22                             <term>theta</term>
23                             <listitem>
24                                 <para>
25                                     [a,b,c] is a real vector of order <literal>3*n</literal>
26                                 </para>
27                                 <variablelist>
28                                     <varlistentry>
29                                         <term>a,b,c</term>
30                                         <listitem>
31                                             <para>
32                                                 <literal>a=[a(1),...,a(n)], b=[b(1),...,b(n)], c=[c(1),...,c(n)]</literal>
33                                             </para>
34                                         </listitem>
35                                     </varlistentry>
36                                 </variablelist>
37                             </listitem>
38                         </varlistentry>
39                         <varlistentry>
40                             <term>p</term>
41                             <listitem>
42                                 <para>(3*n x 3*n) real matrix.</para>
43                             </listitem>
44                         </varlistentry>
45                         <varlistentry>
46                             <term>phi,psi,l</term>
47                             <listitem>
48                                 <para>
49                                     real vector of dimension <literal>3*n</literal>
50                                 </para>
51                             </listitem>
52                         </varlistentry>
53                     </variablelist>
54                     <para>
55                         Applicable values for the first call:
56                     </para>
57                     <programlisting role=""><![CDATA[ 
58 theta=phi=psi=l=0*ones(1,3*n). p=eye(3*n,3*n)
59  ]]></programlisting>
60                 </listitem>
61             </varlistentry>
62             <varlistentry>
63                 <term>u0</term>
64                 <listitem>
65                     <para>
66                         real vector of inputs (arbitrary size). (<literal>u0($)</literal> is not used by rpem)
67                     </para>
68                 </listitem>
69             </varlistentry>
70             <varlistentry>
71                 <term>y0</term>
72                 <listitem>
73                     <para>
74                         vector of outputs (same dimension as <literal>u0</literal>). (<literal>y0(1)</literal> is not used by rpem).
75                     </para>
76                     <para>
77                         If the time domain is <literal>(t0,t0+k-1)</literal> the <literal>u0</literal> vector contains the inputs
78                     </para>
79                     <para>
80                         <literal>u(t0),u(t0+1),..,u(t0+k-1)</literal> and <literal>y0</literal> the outputs
81                     </para>
82                     <para>
83                         <literal>y(t0),y(t0+1),..,y(t0+k-1)</literal>
84                     </para>
85                 </listitem>
86             </varlistentry>
87         </variablelist>
88     </refsection>
89     <refsection>
90         <title>Optional arguments</title>
91         <variablelist>
92             <varlistentry>
93                 <term>lambda</term>
94                 <listitem>
95                     <para>optional argument (forgetting constant) choosed close to 1 as convergence occur:</para>
96                     <para>
97                         <literal>lambda=[lambda0,alfa,beta]</literal> evolves according to :
98                     </para>
99                     <programlisting role=""><![CDATA[ 
100 lambda(t)=alfa*lambda(t-1)+beta 
101  ]]></programlisting>
102                     <para>
103                         with <literal>lambda(0)=lambda0</literal>
104                     </para>
105                 </listitem>
106             </varlistentry>
107             <varlistentry>
108                 <term>k</term>
109                 <listitem>
110                     <para>contraction factor to be chosen close to 1 as convergence occurs.</para>
111                     <para>
112                         <literal>k=[k0,mu,nu]</literal> evolves according to:
113                     </para>
114                     <programlisting role=""><![CDATA[ 
115 k(t)=mu*k(t-1)+nu 
116  ]]></programlisting>
117                     <para>
118                         with <literal>k(0)=k0</literal>.
119                     </para>
120                 </listitem>
121             </varlistentry>
122             <varlistentry>
123                 <term>c</term>
124                 <listitem>
125                     <para>
126                         Large argument.(<literal>c=1000</literal> is the default value).
127                     </para>
128                 </listitem>
129             </varlistentry>
130         </variablelist>
131     </refsection>
132     <refsection>
133         <title>Outputs:</title>
134         <variablelist>
135             <varlistentry>
136                 <term>w1</term>
137                 <listitem>
138                     <para>
139                         Update for <literal>w0</literal>.
140                     </para>
141                 </listitem>
142             </varlistentry>
143             <varlistentry>
144                 <term>v</term>
145                 <listitem>
146                     <para>
147                         sum of squared prediction errors on <literal>u0, y0</literal>.(optional).
148                     </para>
149                     <para>
150                         In particular <literal>w1(1)</literal> is the new
151                         estimate of <literal>theta</literal>. If a new
152                         sample <literal>u1, y1</literal> is available the update is
153                         obtained by:
154                     </para>
155                     <para>
156                         <literal>[w2,[v]]=rpem(w1,u1,y1,[lambda,[k,[c]]])</literal>. Arbitrary
157                         large series can thus be treated.
158                     </para>
159                 </listitem>
160             </varlistentry>
161         </variablelist>
162     </refsection>
163     <refsection>
164         <title>Description</title>
165         <para>
166             Recursive estimation of arguments in an ARMAX model.
167             Uses Ljung-Soderstrom recursive prediction error method.
168             Model considered is the following:
169         </para>
170         <programlisting role=""><![CDATA[ 
171 y(t)+a(1)*y(t-1)+...+a(n)*y(t-n)=
172 b(1)*u(t-1)+...+b(n)*u(t-n)+e(t)+c(1)*e(t-1)+...+c(n)*e(t-n)
173  ]]></programlisting>
174         <para>
175         </para>
176         <para>
177             The effect of this command is to update the estimation of
178             unknown argument <literal>theta=[a,b,c]</literal> with
179         </para>
180         <para>
181             <literal>a=[a(1),...,a(n)], b=[b(1),...,b(n)], c=[c(1),...,c(n)]</literal>.
182         </para>
183     </refsection>
184     <refsection>
185         <title>Example</title>
186         <programlisting role="Example"><![CDATA[
187 nbPoints = 500; // Number of points computed
188
189 // Real parameters a,b,c: here, y=u
190 a=cat(2,1,zeros(1,nbPoints - 1));
191 b=cat(2,1,zeros(1,nbPoints - 1));
192 c=zeros(1,nbPoints);
193
194 // Generate input signal
195 t=linspace(0,50,1000);
196 w=%pi/4;
197 u=cos(w*t);
198
199 // Generate output signal
200 Arma=armac(a,b,c,1,1,0);
201 y=arsimul(Arma,u);
202
203 f1=figure("figure_name","figure1","backgroundColor",[1 1 1]);
204 subplot(3,1,1);
205 plot(t, u, "b+");
206 xtitle("Input");
207 subplot(3,1,2);
208 plot(t, y);
209
210 // Arguments of rpem
211 phi=zeros(1,nbPoints*3);
212 psi=zeros(1,nbPoints*3);
213 l=zeros(1,nbPoints*3);
214 p=1*eye(nbPoints*3,nbPoints*3);
215 theta=[0*a 0*b 0*c];
216 w0=list(theta,p,l,phi,psi);
217 [w0, v]=rpem(w0,u,y);
218
219 // Get estimated parameters:
220 a_est=w0(1)(1);
221 b_est=w0(1)(nbPoints + 1);
222 c_est=w0(1)(2 * nbPoints + 1);
223 for i=2:nbPoints
224     a_est=cat(2,a_est,w0(1)(i));
225     b_est=cat(2,b_est,w0(1)(i+nbPoints));
226     c_est=cat(2,c_est,w0(1)(i+2*nbPoints));
227 end
228
229 // Generate and plot output estimated
230 Arma_est=armac(a_est,b_est,c_est,1,1,0);
231 y_est=arsimul(Arma_est,u);
232 plot(t, y_est,"color","red");
233 xtitle("Real output(blue), Estimated output (red)");
234
235 // Plot estimated parameters
236 subplot(3,1,3);
237 xtitle("a,b,c estimated");
238 plot(a_est(1,:),"color","red");
239 plot(b_est(1,:),"color","green");
240 plot(c_est(1,:),"color","blue");
241 ]]>
242         </programlisting>
243         <scilab:image>
244 nbPoints = 500; // Number of points computed
245
246 // Real parameters a,b,c: here, y=u
247 a=cat(2,1,zeros(1,nbPoints - 1));
248 b=cat(2,1,zeros(1,nbPoints - 1));
249 c=zeros(1,nbPoints);
250
251 // Generate input signal
252 t=linspace(0,50,1000);
253 w=%pi/4;
254 u=cos(w*t);
255
256 // Generate output signal
257 Arma=armac(a,b,c,1,1,0);
258 y=arsimul(Arma,u);
259
260 f1=figure("figure_name","figure1","backgroundColor",[1 1 1]);
261 subplot(3,1,1);
262 plot(t, u, "b+");
263 xtitle("Input");
264 subplot(3,1,2);
265 plot(t, y);
266
267 // Arguments of rpem
268 phi=zeros(1,nbPoints*3);
269 psi=zeros(1,nbPoints*3);
270 l=zeros(1,nbPoints*3);
271 p=1*eye(nbPoints*3,nbPoints*3);
272 theta=[0*a 0*b 0*c];
273 w0=list(theta,p,l,phi,psi);
274 [w0, v]=rpem(w0,u,y);
275
276 // Get estimated parameters:
277 a_est=w0(1)(1);
278 b_est=w0(1)(nbPoints + 1);
279 c_est=w0(1)(2 * nbPoints + 1);
280 for i=2:nbPoints
281     a_est=cat(2,a_est,w0(1)(i));
282     b_est=cat(2,b_est,w0(1)(i+nbPoints));
283     c_est=cat(2,c_est,w0(1)(i+2*nbPoints));
284 end
285
286 // Generate and plot output estimated
287 Arma_est=armac(a_est,b_est,c_est,1,1,0);
288 y_est=arsimul(Arma_est,u);
289 plot(t, y_est,"color","red");
290 xtitle("Real output(blue), Estimated output (red)");
291
292 // Plot estimated parameters
293 subplot(3,1,3);
294 xtitle("a,b,c estimated");
295 plot(a_est(1,:),"color","red");
296 plot(b_est(1,:),"color","green");
297 plot(c_est(1,:),"color","blue");
298         </scilab:image>
299     </refsection>
300 </refentry>