Algoritme voor het schatten van de parameters in ARMA-modellen met meetfouten op de in- en uitvoer

Citation for published version (APA):

Vregelaar, ten, J. M. (1985). Algoritme voor het schatten van de parameters in ARMA-modellen met meetfouten op de in- en uitvoer. (Memorandum COSOR; Vol. 8510). Technische Hogeschool Eindhoven.

Document status and date: Gepubliceerd: 01/01/1985

Document Version:

Uitgevers PDF, ook bekend als Version of Record

Please check the document version of this publication:

• A submitted manuscript is the version of the article upon submission and before peer-review. There can be important differences between the submitted version and the official published version of record. People interested in the research are advised to contact the author for the final version of the publication, or visit the DOI to the publisher's website.

• The final author version and the galley proof are versions of the publication after peer review.

• The final published version features the final layout of the paper including the volume, issue and page numbers.

Link to publication

General rights

Copyright and moral rights for the publications made accessible in the public portal are retained by the authors and/or other copyright owners and it is a condition of accessing publications that users recognise and abide by the legal requirements associated with these rights. • Users may download and print one copy of any publication from the public portal for the purpose of private study or research. • You may not further distribute the material or use it for any profit-making activity or commercial gain

• You may freely distribute the URL identifying the publication in the public portal.

If the publication is distributed under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license above, please follow below link for the End User Agreement:

www.tue.nl/taverne Take down policy

If you believe that this document breaches copyright please contact us at: openaccess@tue.nl

providing details and we will investigate your claim.

Memorandum-COSOR 85-10

Algoritme voor het schatten van de

parameters in Arma-modellen met

meetfouten op de in- en uitvoer

door

J.M. ten Vregelaar

Eindhoven, Nederland

O. Inhoud 1. 2. 3. 4. 5. 6. 7. 8. 9. Inle1ding Probleemstelling

Bet schatten van de parameters

Algori tme m. b. v. Q - R decompositie·

Q - R decompositie gebaseerd op Householder transformaties

Bepaling van A uit RA =u

Bepaling van

n ,

t

af

Bepaling van de gradient

aYk Conclusies 10. Referenties 1 2 5 9 12 23 25 27 30 33

1. Inleldlng

In dlt memorandum wordt een algoritme beschreven voor het schatten

van de parameters in Arma-modellen met meetfouten op in- en uitvoer.

We be schouwen hier het MIMO (multiple input - multiple output) geval

en nemen aan dat het totaal aantal waarnemingen van de in- en

uitgan-gen (= N) groot is vergeleken met het totaal aantal parameters

2

ps + (q + 1)sr.

In [EIS] staat een algoritme beschreven voor het SISO-geval (single

input - single output), dat gebruik maakt van zgn. Toeplitz-afstanden

om een grootschallge matrix te inverteren. In de hier gehanteerde, meer

inzichtelijke methode, komt een Q-R decompositle in de plaats van

ma-trixinversie.

Na de probleemstelllng in paragraaf 2 voIgt in de paragrafen 3 en 4 de

methode die aanleiding geeft tot het algoritme. In het algoritme

onder-scheiden we 4 stappen die achtereenvolgens in de paragrafen 5

tim

8 aan de orde komen. De cruciale stap is de Q-R decompositie, waarvan aan

het eind van paragraaf 5 een algoritme in pseudo-algol vorm is opgenomen.

Arma-modellen met ruis op de in- en uitgangen vinden o.a. toepassing in

2. Probleemstelling

~t afmeting r

'"I

G c.k afmeting·s

We beschouwen een dynamisch lineair tijdinvariant MIMO (multiple input

-multiple output) systeem, dat als voIgt gemodelleerd wordt (zgn.

ARMA-model) :

+ ... + +

Auto Regressief deel Moving Average deel

t

=

1,2, .••• De dieptes p en q zijn bekend (P.q ~ 0).

De parameters a

1 en S1 zijn s x s resp. s x r matrices.

Er z1jn met fout gemeten waarnemingen beschikbaar:

x

=

~ + et

t t t t = 1,2 ... N .

We maken de volgende aannames voor de fouten e

t en e~: lEet _{= 0} t var et t (1) (2) (3)

op grond van de waarnemingen (2) waarbij de ware waarden tt' ~ van

de 1n- en uitgangen voldoen aan het model (1).

We kunnen model (1) ook weergeven ala

A'I1 + B~ = 0 , (4) waarbij -I

₁

_So

•

_•

_•

:ll

•

₁₁₁ 0. 1

•

• •

g

•

•

•

_•

g

•

•

_•

• •

_•

_•

•

•

•

•

,

A

=

•

_•

B

=

•

_•

11

=

~

=

~J

•

•

,

,

•

_•

•

_•

_TN

•

-

•

-

•

•

a

-

•

_•

_Sq.

•

p.

_•

•

_•

_•

•

_•

•

•

_•

_•

•

•

_•

_•

g _a p

...

~

-I g

-13 •••

q

So

sN x sN sN x rN

A en B zijn band, llnksonder blok Toeplltz matrices. Bovendlen is A

regulier.

Na het op analoge wijze invoeren van de kolomvectoren y,x,e en e'

no-teren we (2) en (3) als y = '11 + e x

=

~ + et _(5a) lEe

=

0 Ee'

=

0 (5b) e,e' onafhankelijk (5c) en

reJ

2 VAR

b'

=

(j I (5d)

Opmerkingen.

1. We veronderstellen dat voor grote N het begineffect te

verwaar-lozen is, waardoor (4) gerechtvaardigd Is. Eventueel kunnen ;t en n

t voor t ~ 0 geschat worden.

2. Aanname (5d) garandeert de identiflceerbaarheid van de parameters,

3. Het schatten van de parameters

We bepalen schatters voor ai' Sj' n en ; door de restkwadratensom te

minimaliseren:

Ilx -

;112

+ 111 - nl1

2

onder An +

Bt

= 0 . (6)

Als we aannemen d.t de fouten e en et _{normaal verdeeld zijn dan zijn de}

oplossingen.van (6) de maximum likelihood schatters voor de parameters.

Na het invoeren van de notaties

gaat (6) over in

min liz -

011

2 onder Do = 0 • 'V,eS

Omdat A regulier is, is D van volle rijenrang.

Voor vaste 'V geeft toepassing van Lagrange

1 T T

L(O,A)

=

2(Z - 0) (z - 0) + A Do

we lossen 0 en A op uit VL

=

0, dus uit

'{o

+

Do = 0

(7)

(Sa)

Dan is

Dus z - 15

=

D Dz +

Probleem (7) is nu gereduceerd tot

min fey) met fey) = IID+DzI12

y

(9a)

+

( l - D D)z (9b)

(10)

Tot slot van deze paragraaf leiden we uitdrukkingen af voor de eerste en

tweede afgeleiden van fey), die we eventueel benutten b1j de oplossing van

probleem (10).

()f

Noteer fk := aYk met Yk is een element van een parametermatr1x, analoog

Dk , ).k' ok' Vi t (8) voIgt Ok + DT).. _k

=

_DT ). _k (l1a) DOk

=

-D 6 _k (llb) Los Ok en Ak op uit (11) : + DDT). (l1a) _ DDT ). DDT). (Ub) DDT A DOk _k

=

_k dus _k

=

_{D k}

° -

_k

waaruit voIgt

(12a)

(12b)

Eerste afgeleide

Differentieer f(y)'

=

~z

_

a~2:

(13)

voor Yk afkomstig uit een ai-~trix is fk

en voor Y_k uit een Sj-matrix is

Tweede afgeleide Noteer dan is • Dk constant

=

AtTDka + ATDka t (1·

Opmerking

1. In termen van A, Bt Y en x is £(y) te schrijven a1s (zie (10»,

T

T

T -1

£(y)

=

(Ay + Bx) (AA + BB) (Ay + Bx) • (15)

2. Uit (8a) voIgt dat 0 naast (9b) ook geschreven kan worden ala

dus 11

=

Y

(16)

4. Algoritme met behulp van Q - R decompositie

We beschrijven in deze paragraaf een algoritme dat bij gegeven y de

waarde f(y) en de gradient

~~

berekent. Bij geschikt gekozen start-waarden, zal een numeriek minimaliseringsalgoritme dat gebruik maakt

van de gradient het gezochte minimum bepalen.

We hebben problee~ (6) herschreven tot:

en

verder is

min f(y) , met f(y)

=

Ilo+oz

112

y

als Yk afkomstig uit een ai-matrix

af

aYk

=

als Yk afkomstig uit een Bj-matrix ,

).

=

(O+)Tz

n

₌

y - AT).

E;

=

x _ BT).

In [EIS] wordt voor elke waarde van y voor het SIS~geval een matrix geinverteerd m.b.v. zgn. Toeplitz-afstanden. Een alternatief wordt

ge-geven door [GOL], p. 410 e. v.: Q - R decomposi tie. In de volgende

para-(17a) (17b) (17c) (17d) (17e) (17f)

graaf komt de uitvoering van de Q - R decompositie uitvoerig aan de orde.

Hier gaan we in op de consequenties ervan voor de berekening van f(y)

at

en

-a- .

Gelet op de nevenvoorwaarde Do = 0 (zie (7)) maken we een Q - R

decompo-. T

sitie van de matrix D die van volle kolommenrang is:

DT

=

Q[~]

, Q is orthogonaal, R 1s rechtsboven, band en regulier. Als Q

=

(Ql

I

Q2) voIgt T D

=

Q 1R . Dan is en

Uitgedrukt in z,

Q,

R voIgt voor f en

A

(zie (17.a,d»

Definieer dan 1s en T u := Q 1 z f(y) =

IIul1

2 -1

A

=

R

u • (18) (19a) (19b) (20a) (20b) (20c) •

_af

T (18) T T

V~~r de berekening van f(y) en

-a--

hebben we nodig u

=

Q₁ z, R = Q

1 D

Yk -1

= R u, waarbij de waarnemingsvector z gegeven is.

en

A

Breid de matrix DT uit met de kolomvector z dan geldt

(21)

Aan de hand van bovenstaande afleiding komen we tothet volgende

algo-af

ritme voor de berekening van fey) en --- (y) als y gegeven 1s. aYk

1. Maak een Q-R deeompositie van DT; via (21) volgen R, u en f

=

Ilu112,

2. Los X op uit RX

=

u, R is reehtsboven,'band ..

3, Bereken n =y-.AX T en 4, Bereken T t

=

x -

B A

(zie (17 e,f». (de (17 b ,e» .

In de volgende 4 paragrafen w~rden de 4 stappen in het algoritme

5. Q - R decompositie gebaseerd op Householdertransformaties

5.1. Algemeen

Zij A een m x n matrix (m ~ n). van volle kolommenrang.

Dan 1s er een orthogonale matrix Q (zie b.v. [GOLJ, p. 147 e.v.) zodanig dat

R is rechtsboven en regu11er.

Door de eis dat de diagonaalelementen van R pos1 tief zijn is deze Q - R

decompositie van A eendu1dig bepaald.

Een methode ter berekening van Q en R wordt gegeven door gebru1k te maken van Householder matrices (zie [GOLJ, p. 38-43).

n T T

Als v E lR , v rF 0 dan is P : = I - 2vv Iv v een n x n Householder matrix. Het is duidelijk dat P symmetrisch en orthogonaal is. De voor ons doel

beoogde toepassing van Householder matrices is de volgende.

Stel x

e:

mn

(22)

Opmerking. Om cijferverl1es tegen te gaan, nemen we v = x + sign(x

1) IIxllel (zie [STEJ. p. 233).

Onderstaand voorbeeld laat zien hoe we Q en R berekenen voor gegeven A.

Voorbeeld. Stel m

=

6, n

=

5 en Householdermatrices P

1 en P2 zijn berekend zodanig dat

x x x x x 0 x x x x P 2P1A 0 0 x x x

=

0 0 x x x 0 0 x x x 0 0 x x x

Bepaal nu een 4 x 4 Householder matrix ~3 zodanig dat

x x x 0 ,...

_...

P 3 x

=

0 en neem P3

=

diag(I2 ' P3) x 0 dan is x x x x x 0 x x x x 0 0 x x x P 3P2P1A = 0 0 0 x x 0 0 0 x x 0 0 0 x x Pi is symmetrisch en orthogonaal, i

=

l, ... ,n. R

'5.2. Speciale Q - R decompos it ie

T

We willen nu bovenstaand recept toepassen op de matrix D , Om de Toeplitz

T

structuur in de matrix D niet kapot te maken, voerenwe een speciaal op

die structuur toegesneden Q - R decompositie uit.

Volgens (21) zijn we geinteresseerd in de werking van Q op (D Iz). T T

T .

We schrijven Q als produkt van Ns matrices, die alle orthogonaal en

symmetrisch zijn:

QT

=

H

1 N,s

H H H H

1 - 1

2_,_s _ _ _

2_,_11 1..-,1_

1 ,_s---,. __

1_,_11

(23)

We zullen nu d.m.v. de werking van H

1•1 laten zien hoe de hlok Toeplitz-structuur behouden b11jft.

De werking van Hj,t in het algemeen (j

=

l, ... ,N; t

=

l, ... ,s) komt daar-na aan de orde. Definieer T a₁ := a i b i Sf T := m := max(p,q) i

₌

i

₌

1, .•. ,p 0, •.•• q .- -I s T

dan kunnen we, door eventueel enkele blokken a₁ of b

i nul te nemen, (D Iz)

a

_{o ...}

a Y1 • m· • . • •• a ••

.

m

•

_•

• •

Vermenigvuldiging met H 1,1 '11

Beschouw de eerste kolom in ( D

I

z). Deze bevat slechts (r + 1) elementen die ongelijk aan 0 zijn: a

O ' b

o

, ...

,bO

1,1 1,1 r,l

....,

Laat h een (r + 1) x (r + 1) Householder matrix zijn, waarvoor geldt dat

a

_*

0 1 1

,

b 0 ...., 0 1 1 h

,

=

We kunnen nu dezelfde hook op de overige (N - 1) onderelkaarstaande

a

_o

en b

h11 h12 • •

•

_•

h11 H1 ,1

=

_T h12 h22

•

_•

_•

• T h12 waarbij

t!!1=

diag(~1,1·1

•.••• 1) SXS

~

fi2

_•

2

·

.

·

.

·

.

,...,

,...

h2 ,r+l···hr+1,r+1

Ii'

_2,r+1

•

_•

h12

•

_•

h22 .::...l (s + r)N x (s + r)N ,

Het gevolg is dat in H 1,l(D

T

lz) aIle eerste kolommen van de blokken op

de bo-diagonaal nul geworden zijn, en dat de overige (8 - 1) kolommen

van die blokken identiek getransformeerd zijn. Zo transformeren aIle

onderelkaarstaande paren a

1, b1 tim am' bm op identieke wijze, zodat

T

de blok Toeplitz structuur in de matrix D behouden blijft.

We kunnen voor aIle volgende vermenigvuldigingen een analoge redenering

ophangen, zie "vermenigvuldiging met H

j f t , Essentieel is dat de matrix BT

T

en het "resterende dee1 tt van de matrix A blok Toepli tz blijven door de

speciale keuze van de Hj,l-matrices.

Vermen1gvuldiging met Hj,t

We hoeven slechts de werking van een (r + 1) x (r + 1) Householder

matrix

h

op onderstaande blokken en vectoren te berekenen om de ttnieuwe" (DTI z) helemaal teo kennen. Door het behoud van de Toeplitz structuur

blijven de blokdiagonalen constant (in de matrix AT aIleen de je tIm Ne

blokrij). De veranderde elementen behouden hun "naam":

a 1 b. l-b i en a a i it t t,s

,

l,s b b i i 1 ,t 1,s _{:= h}

....

l,s i=O, •.. ,min(m, b i b. r,t 1-r,s r,s 1s2 2 .2

I

-sign(a O ) aO + bO + •.• + bO t,t t,t 1,t r,t = b O := 0 (zie (22,25» r,t en

~:::'j

:=

h

~::'j

k = O,l ••.. ,N - j .

...

Hierin is h Householder, .... T T h

=

I - 2vv Iv v , met N - j) (24a) (24b) (25) + •.• + b2

_o

I

_'

b _{0 , . . . , 0} b ) T r,t 1,t r,t

De volgende keuze voor de matrices H

j •t in (23) zorgt voor handhaving

van de Toeplitz structuur in de matrix DT:

[1

~ (j-l) I •

~I

h22 j = I •...• N; t=1 (s+r)Nx(s+r)N waarbij

ful._

..., diag(I, ...• 1, hI 1,1, ..•• 1) SXS ,-

_,

t e 1 plaats "'" _h 2,r+l ,.., t v h2 _,r+ 1" h _{r+ ,r+} 1 1 • • h22 I

•

~I=

_sXr

•

_•

r

I

[1.2

~ ,.., .•• hI _.r+ 1 ~ e De werking van H

j•t bestaat uit het schoonvegen van de 1 -kolom in aIle

-+- tE

b_O blokken op de blokdiagonaal N+l,j; N+2, j+l;, •• ; 2N-j+I,N. Naelke H

j

'I is er een blokdiagonaal met bO-blokken nul geworden,

T T

Schema: Werking van Q op (D Iz)

Het volgende schema beschrijft globaal de werking van QT op (DTlz); we

definieren daartoe R := (r

ij), r1j 1s (s x s)-blok voor 1,j = 1, ... ,N

I I

1 en U := met u i (s x l)-vector. uN r-

-r '" , .... 1,1 r 1 , m+1 u1 ~ .••.. 0 000 a _{m Y2} •

•

_'

••

...

..

a m

..

_•

• •

_..

•

P1(DTlz) • a _Y N 0 =

~

0 0 0.0 t • • • " 0 b _xl m-1

.

_•

•

_b m-l

•

_•

•

_•

• •

'b 0 ~-l 0 w N

-r-

-r 1 1

•

.

_..

o .• or 1 ,m+

_•

1 u1

..

•••

_..

.

..

_'

.

...

_u r _{N-m, -m} N0' 0 or _{N-m N N-m}

~

_,

'.0 a m-1 YN- m+1

·

•

_..

_·

0

·

0

•

_..

0

_·

a

_o

_Y~

~

.. o b _m-1 xl

,

....

,

•

_•

_·

0

•

0

_·

•

_•

_·

0 b

_o

x m W m+1

·

·

·

w N

-

-r

1 1 ...

,

R en u zijn nu bepaald, (zie 21).

Toelichting bij het schema

•

_•

•

_•

r 1,m+1

• •

'

.

• •

_•

• •

_·r N-m,N

•

(26) 1. We hernummeren de bi-blokken is na elke vermenigvuldiging

na elke vermenigvuldiging met H

j ,s . Verder

e e

met H

j ,s de j blokrij van R en de j

blok-component van u gevonden:

rj,j := aO····, rj,j+min(m. N-j)

.-en

Deze b1ijven onveranderd. evenals w l '

N-j+

a

min(m, N-j)

e

2. Omdat voor j

=

2, ... ,N-m het aantal blokken a

i in de j blokrij een meer

e

bedraagt dan het aantal blokken b

i in de (N + 1) blokrij van (de

getrans-T

formeerde) D blijft het aantal blokken b

i constant (= m) voor j = 2, ... ,N - m enerzijds verdwijnt het b

_o

blok, anderzijds wordt het bm-blok gevormd

dat na hernummering de naam b 1 krijgtj in (24 a) is we11swaar het b

-m- m

blok 1n het rechter11d nul, maar het am-blok niet.

3. Omdat na elke stap a

O rechtsboven bl1jft 1s slechts een (r + 1) x (r + 1) ondermatrix van de (s + r) x (s + r) mat:1>;. die op

~:J

werkt van belang.

De (r + 1) x (r + 1) Householder matrix h (in (24» wordt aangevuld met

s - 1 termen 1 op de diagonaal .

In bet bizonder zijn aIle r

1

_•

1 (1

=

1, ... ,N) in (26) recbtsboven . "'" T T 4. U1 t h = I - 2 vv Iv v voIgt

h.

j

=

6_i . 1, , J 2 - - - v v Ilv112 ' 1 j 1,j = 1, . . . ,r+l •

Bovenstaande bescbrijving geeft aanleiding tot bet volgende algoritme ter

bepaling van R,~ en fey).

Invoer aO, .. ·,a_m, _{bO,···,bm, z} =

[~J

Ultvoer R, u, f •

Bescbouw bet geval m ~ 1. Het commentaar is vermeld tussen accolades. De

waarden voor j en t corresponderen met de vermenigvuldiging met H

begin f := 0 ;

end

for j := 1(1) N do

begin if 2 :5 j S N - m then b :

=

0

m {zie toelichting punt 2}

end for 1 := 1(1) S do begin v := (a O ' bO i,l I fi {de (25)} {~ie (24 a)} for i := 0(1)m1n(m,N - j) do "verander a 1 ,bi" 1*

for k := O(l)N-j do "verander Y

k+j ,t;Xk+l" {zie (24b: for i := 0(1)m1n(m,N - j) - 1 do begin b i

:=

bi+1 {Hernummer de bi's} r j. j+i := a₁ r := a

j,j+min(m,N-j) min(m,N-j) {zie toelichtil

6. Bepaling van A uit RA

=

u

Via (2.6) schrijven we RA

=

u als

r:

_{'&'1,1 ... l,m+1} r

.

'.

_'. r N-m,N A1 u_l = AN

lJ

met r₁,j : s x s , Ai,u_i : s X l . e

Noteer Ai,j en Ui,j voor de j component van de vectoren A1 resp. U_i

(1 = 1 •••.• N; j = 1 •..• ,s), z1e ook par. 7 en 8.

R en U z1jn bekend (zie paragraaf 5), R is een rechtsboven en

band-matrix, ook alle r

i

,

i zijn rechtsboven voor 1

=

1, ... ,N.

Schrijf (27) uit door de N blokrijen te vermenlgvuldigen met A, te

be-ginnen met de Ne:

rN,NAN

=

UN

r N- 1 ,N-IAN-1 = _{uN_1 - rN_1,NAN}

r A =u - r A -N-m,N-m N-m N-m N-m,N-m+l N-m+l r A N-m-l,N-m-l N-m-l = u._ N-m-l . - r N-m-1'N-m N-m A - - r N-m-1,N-l N-1 A - r A • l.m+l m+l (27) e Los uit de eerste vergelijk1ng AN op.uit de tweede A_N__l " .. , uit de N Al'

Algoritme.

min(N,i+m)

__ I.2!:

i:=N(-l) 1 do begin rechterlid :=u

i -

1:

r1 kAk

k=i+l '

rechter1id"

Het dee1probleem los A op uit rA = u, 1s hetze1fde probleem met dit

verschil dat r nu een s x s matrix is en

A

en u s x

1

vectoren, en

verder 1s r weI rechtsboven, maar geen bandmatrix.

uitgeschreven geeft r

A

s,s s

=

Us r

A

=

u - r

A

s-l,s-l s-l s-l s-l,s s - r A-I,s s

Het bijbehorende algoritme luidt nu:

s

for i := s(-l)1 do begin rechterlid := u_i -

I

r 1kAk k=1+1

rechterlid / r

i

,

1

7. Bepallng fl, ~

Analoog aan A _{en u noteren we 1'11 ,j' l;i,j' Xi,j en Yi,j voor de je} de e ( i

=

1, ... ,N en j = 1, ... , r s) .

component van i vector resp

We berekenen fl en l; uit (zie (17 e,f»

(1'1 = Y - AT).. (28a)

i.

= x -

i?A

(28b)

V~~r

elke y _{= (al,···,ap '} 6_{0 ,···,13q ), ligt}

[~J

vast; z

=

[~J

is gegevenj voor A en S zie (4) en A is het resultaat van de berekenlnge~ uit

para-graaf 6.

1-1

• a₁ T a T Al ...

.

p, T a AT" .p T min(l,N-i) T =

.

,

_{dus fA A)ll= -A +} a A . , i k=l k k+1 , _T " _, a₁ sxl , _, " _{-I AN i}

₌

_{1, ... ,N}

.

Iso

T ," .

a

_{q ...} T Al

...

T

'a

sTA

₌

.

.

.

.

q _dus T min(I,N-i) T

I(S A)

11

= 13 A

"

_"

_...

,

k=O k k+i

,

. e

T r xl

0 _AN

U1t (28) voIgt nu 1

=

1, •..• N; j=l, ...• s

~1,j

1 = 1, ...• N; j=1, ...• r Omdat T S T s (a_{k Ak+1 ) =}

I

(a_{k )} j 1 A

=

I

(a k) hk+1 ,1 j = 1, ...• s j 1=1 • k+i,1 1=1 1,j en T s s (Sk A_{k+1 )}

=

I

(a

T) A =

I

(a

_k) _{Ak+1,1 j}

-

I, ... ,r j 1=1 k j,t k+i,1 t=l 1, j v1nden we voor n en ~ (29a)

n

1 _.J . min(i!,N-i) s

=

Y1.j + A₁,j -

L.

I

( a ) A i 1

=

1, .•• ,N: j=1,., . k=l 1=1 k 1. j k+ ,1 (29b)

~i.j

i

=

1 •... ,N; j=l,.

8 . B epa ng van e gra ent -a--Ii d di " af Yk

We berekenen de gradient ; f uit de formuies (17 b,e)

Yk

af aYk =

als Y_k afkomstig u1t een ai-matrix

alsY

k afkomstig uit een 6j-matrix •

De veetoren A, n en ~ zijn bekend (u1t paragraaf 6 en 7), en

(30a)

(30b)

We voeren eerst de vermenigvuldiging

AT~

resp. ATBk uit, dit levert sN-resp. rN-rijveetoren.

T

Bedenk dat A

Na 3 speciale Yk'S voIgt het algemene resultaat.

Dan 1s 0.8.+1 : [ Bk

=

1.s + i ~ (N':l)s + i -+ 1 t O.r+j 1.r+j •.• (N-1)r+j

:

--.-~

]

Bk be vat N elemente~ 1, de overlge zijn O. Gevolg plaats: t j t r+j t (N-1)r+j 1 - 1, •..• s j = l , .•.• r .

De over1ge elementen in deze rN-rijvector zijn O ..

2) Y_k = (a 1) 1.j s+i -+

[:

]

Ak = 2s+i

.

.

-+ 1 ••

...

(N-l)s+i -+ ~ · · . 1 t _t j s+j (N-2)s+j Ak bevat (N-1) elementen 1. Nu 1s 0) 1, j :: 1, ...• s . plaats: t j t s+j 3) Y k

=

(61) analoog aan 2). 1,j "TB

₌

_{("2 i 1..3 1} k

,

_•

t t plaats j r+j Algemeen geldt: • t (N-2)s+j AN 1 0)

•

t (N-2)r+j 1

₌

1, ..•• s j

=

1, ..•• r

AT

3A = °k+l, i Ak +2 ,1 AN 1 0 0) 3 (a. k)

..

,

1,j t t t plaats j s+j (N-k-l)s+j

AT

aB ("k+l,l "k+2,1 AN i 0 0) =

..

a(Sk)

,

1,j t t t j r+j (N-k-l)r+j

Berekening van de afgelelde uit (30) resulteert in

N

=

2 ~ n=k+l N

l:

n=k+l A

n

n,i n-k,j , A

t

n,i n-k,j . k = l , •.• ,p 1 , j = 1 , ... ,s k=O, ..• ,q 1 = l , •.• ,s j = l , •.• ,r 1,j=1, ... ,s k= 1, ... ,p

.

1=1, ... ,s j = 1; ..• ,r k=Op .. ,q

.

(31a) (31b)

9. Conclusies

We geven een samenvatting van de hier beschreven methode die tot een

algoritme leidt voor het schatten van de parameters in Arma-modellen

met meetfouten op in- en u1tvoer.

We bepalen de schatters voor de parameters y door de

residuenkwadraten-somte minimallseren met de modelvergelijking als nevenvoorwaarde (z is

de waarnemingsv~ctor van de uitgangs- ingangskolom 0):

min liz - 0112 onder Do

=

0, D(y) van volle rijenrang. y,o

Met behulp van de Lagrangiaan L

=

~"z

- 0112 + XTDO v1nden we het probleem zonder nevenvoorwaarde

min fey) met fey)

=

IID+DzI12, D+ = DT(DDT)-1 •

y

Omdat we f niet expliciet als functie van y kennen, lossen we dit

pro-bleem numeriek op. Bij geschikt gekozen beginschattingen zal een numeriek

algori tme dat gebruik maakt van de afgelelden

o

=

z

-ons schattingen van de parameters leveren. Vanwege de vervelende termen

T -1 -1 -T + T

(DD)

=

R R en D D

=

Q

1Q1 in (14) benutten we de tweede afgeleide niet.

T

Als we een Q - R decomposit·le van D geconstrueerd hebben (Q

=

(QIQ2» . dan wordt

f

=

IIul1

2 , metu

en

A

=

R -1 u.

Het algoritme voor de berekening van f en fi als functie van y ziet er

globaal als voIgt uit:

Door een speciale Q - R decomposi tie gebaseerd op Householder matrices

T

wordt de Toeplitz-structuur in de matrix D goed gebruikt.

2.

Los

A

op uit

RA

= u.

R

is een rechtsboven en bandmatrix. Het 1s

moge-lijk de blokstructuur te handhaven bij de terugsubstitut1e. Over1gens

verwachten we dat dit probleem goed geconditioneerd is als de

parameter-waarden 1n DT niet te groot zijn.

3) Bepaal 0 =

z -

D A • T

2

(1 = It''''ps + (q + l)sr) .

Door de bewerkingen op een geschikte wijze te organiseren komen bij de

laatste twee stappen de formules (29) en (31) tevoorsch1jn.

Merk op dat nergens expliciet gebruik is gemaakt van vrije parametrisatie

parametrisaties, die zich in o.a. economische toepassingen voordoen.

Alleen formule (31) voor de afgeleide behoeft enige aanpassing.

Of en zo ja, hoe snel, een computerprogramma gebaseerd op bovenstaand

algoritme ons parameterschattingen levert heeft nu onze interesse. We

zullen dan de resultaten vergelijken met die van het algoritme

be-schreven in [EISJ; eventueel verdient de keuze van de beginschatting

nog wat aandacht.

Het toekomstig onderzoek richt zich verder op de asymptotische

eigen-schappen van de schatters en op toepassingen in de vorm van econometrische

modellen.

In een later stadium zullen we het geval bekijken dat enkele van de

10. Referenties

[EIS]

[GOL]

[STE]

Elsing F., Linssen H.N., Rietbergen H.,

System identlfication from nolsy measurements of inputs and

outputs, System and control letters 2, 348-353, 1983.

Golub G.H., van Loan C.F.,

Matrix Computations, North Oxford Academic, Oxford, 1983

Steweart G. W • ,

Introduction to matrix computations, Academic Press, New York,