Opravy závislých údajů soustavy čidel se známým rozdělením chyb m~řenf

681.586

Ing. Jiří Tůma, CSc. TATRA - kombinát, Kopřivnice

Článek analyzuje p'ostup opravy {uJ,ajt"" čidel pro zadaný zá.kon rozdělení náhodnýoh ohyb měření a pro zadanou sOU8tavu lineárníoh rovnio, které vážou správné hodnoty měřenýoh údaj1i a jejichž počet je menší než počet údajů v jednom pozorovánL. Je od!Jozen obeoný algor'itm1~s výpočtu oprav,. kte,.,,) je konk1'etizován pro normální a rovnomě'r11.é rozděle7!.í ohyb měření. J e d~skutován rovněž postup oprav parametrů hustoty rozdělení ohyb při opakovaném pozorování.

KlíČová . slova: opravy ohyb měření, normální rozdělení ohyb, rovnoměrné rozdělení ohyb, opravy parametrů h1~.stoty rozdělení ohyb měření.

1. Úvod

Věrohodnost údajů soustavy čidel je možné prl splnění. určitých podmrnek posuzovat podle toho, zda vyhovují. matematickým rovnicí.m, které popisují. zákonito'sti sledovaného proceEOl. Mezi tyto zákoni­tosti patří zákon zachováni hmoty a energie, Kon­krétně například součet proudú směřujících do uzlu musí být nulový. Podmínkou výužití této rovnice. ke kontrole měřených údajů je měřeni všech proudů vstupujících do uzlu. Rovnice, kterým mají sprá.vné hodnoty měřených údajů vyhovovat, jsou zpravidla vzhledem k těmto hodnotám lineární nebo je lze snadno linearizovat. Vazebnich rovnic pro měřené údaje múže být několik. Další podmí.nkou ovšem je, aby tyto rovnice byly vzájemně nezávislé a jejich počet byl menší. než počet měřených údajů. Nezávis­lost soustavy vazebnich rovnic je chápána tak, že žádnou rovnici soustavy nelze odvodit z rovnic ostat­ních. Změřené údaje obvykle výše uvedeným rovni­cím nevyhovují, protože jsou zatíženy náhodnými a systematickými chybami měřeni.

Pro opravy závislých pozorováni jsou známy meto­dy z vyrovnávacího počtu . Aplikaci na bilanční vý­počty procesů chemické technologie popisuje Hlavá­ček a kol. [1]. Měřené údaje jsou oprl.tvovány na zá­kladě jednorázově zvolených hodnot parametrů hus-

. toty rozděleni náhodných chyb. Opakované měřeni však poskytuje možnost zpřesnění. statistických cha­rakteristik chyb měřeni. Návrh konkrétnich algoritmů výpočtu vyžaduje podrobný rozbor jejich teoretic­kých základů a kritérií, kterými se mají opravy ří.dit.

2. Definice úlohy

ných údajů vektorem v = (Vl"" Vn)T . Vektor chyb je definován tak, že platí. e = u-v.

Vzájemné vazby mezi správnými hodnotami měře­ných údajů lze vyjádřit soustavou lineárních rovnic

A . v = bo (1)

Matice A je obdélniková typu m X n a vektor bo má počet složek m, přičemž m < n. Po dosazeni za vektor v nabude soustava rovnic (1) tvar

A.e= b (2)

kde b = A. u - bo. Rovnice soustavy (1) mají. být podle zadáni vzájenině nezávislé. To znamená, že ze sloupců matice A lze sestavit regulárni matici řádu m. Hodnost matice A jB h(A) = m. Soustava rovnic (2) vzhledem k většímu počtu neznámých než rovnic nemá jednoznačné řešení.

Chyby měřeni jsou náhodné veličiny. Jejich sta­tistické vlastnosti popiEOlje sdružená hustota rozdě­leni, která je označena p(x), kde x je vektor z n-roz­měrného euklidovského. prostoru. Kromě vektoru x je hustota rozděleni funkcí řady parametrů, daných konkrétním typem rozdělení náhodných chyb. Oprava měřených dat bude v další.m textu odvozena pro nor­mální. a rovnoměrné rozděleni náhodných chyb měření.

Příklad 1

Nechť je dvěma různými způsoby měřen jeden prů­tok. Napřiklad integrovaný údaj ~ kontinuálniho měřiče průtoku je doplněn údajem o sumárním prů­toku U2• Platí

(3) Nechť soubor měřených údajů (jedno pozorováni) je dán vektorem u = (~, ... , u,.)T, chyqy měřeni tj. vektorem e = (Bl' ... , Bn)T a správné hodnoty měře- A = (1, -1), b = (~-U2) .(4)

" ",

-=3-=-1 --,(-=-19:c:S..:..S,-) --,-č_, _10 ___________________________ 259 automatizace

Opravy měřených údajů pro konkrétní tvary sdružené hustoty p(xv x2) budou určeny v příkladech k násle­dujícím kapitolám,

3. Obecné vzorce pro odvozeni algoritmu výpočtu oprav

Vektor chyb s může před uskutečněním měření na­být libovolné hodnoty v souladu se svým rozdělením, Akt měření váže realizaci náhodné veličiny s podmín­kou (2), Volba opravy měřených údajů s může být

libovolná, protože soustava (2) nemá jednoznačné řešení. Jestliže js.ou jednotlivé s.ložky vektoru s

vybrány tak, aby bylo rozumné posuzovat euklidov­skou délku rozdílu s - s*, pak z možných voleb opravy s* mohou být některé vhodnější než jiné [2]. Protože informace o konkrétní hodnotě s je nahrazena sdruženou hustotou rozdělení p(x), lze hodnotit vhod­nost výběru pouze podmíněnou střední hodnotou

J = E{(S-S*)T ,(S-S*) lAs = b} (5)

kde E {.} je symbol operace střední hodnoty, Toto kritérium pro výpočet nejvhodnější opravy s* je možné upravit do tvaru

J = E{(s - me)T. (s - mel I As = b} + + (me-S*)T.(me-S*) (6)

kde me je podmíněná střední hodnota vektoru chyb měření

me = E{s I As = b} (7)

Prvý člen pravé strany rovnice (6) je rozptyl přípust­ných hodnot chyb měření s od střední hodnoty m,. Druhý člen pravé strany (6) je druhá mocnina eukli­dovské délky vektoru m. - s*. Minimální velikosti nabývá kritérium (6) pro s* = m.. Nejvhodnější opravou vektoru pozorovaných údajů podle kritéria (6) je tedy vektor podmíněných středních hodnot m •.

Pro odvození vzorců k výpočtu vektoru m. bude v systému rovnic (2) rozdělena matice A na bloky a vektor s na složky následujícím způsobem

[A l' A 2] [::] = b

kde A l je matice typu m X (n - m) a A 2 je čtvercová matice řádu m . Pořadí složek vektoru s musí být voleno tak, aby matice A 2 byla regulární, tj. hodnost h(A 2) = m a existovala A2 l. Smyslem úprav je zave­dení transformace vektoru s na vektor '1], která zjedno­duší zápis vazební podmínky (2) pro chyby měření. Vhodnou transformací je

'1] = Q .s (9) B maticí

(10)

ve které E označuje jednotkovou matici a O nulovou matici. Vektor '1] je rozdělen na složky '1]1 a '1]2' J edno­duchým výpočtem lze ukázat, že podmínka (2) je ekvivalentní

(ll)

Akt měření tedy znamená přiřazení konkrétní hod­noty složce '1]2'

Náhodnému vektoru s , který je l10zdělen na složky Sl a S2' přísluší sdružená hustota rozdělení

(12)

s podobně rozděleným vektorem x na složky Xl a x2 .

Sdružené hustotě rozdělení náhodného vektoru '1] se složkam'i '1]1 a '1]2 přísluší

(13)

Hustotu rozdělení g(y) určuje transformace (9) s ja­kobiánem rovným jednotce

g(Yl> Y2) = P(Yl' Y2 - A21~ Al Yl) (14)

Výpočet oprav dospěl do fáze , ve které je odvozena sdružená hustota náhodného vektoru '1] se složkou '1]2 nabývající po aktu měření známou hodnotu (ll). Je třeba určit, jak znalost velikosti jedné složky ná­hodného vektoru '1] ovlivní sdruženou hustotu rozdě­lení druhé složky, Z axiomatických základů teorie pravděpodobnosti plyne

gl(Yl I '1]2 = A2:1b) = g(Yl' Y~) I -1 (15) g2(Y2) Y2=A. b

kde (16)

je marginální sdružená hustota rozdělení složky '1]2' která se získá integrací (13) přes obor hodnot proměn­né Yl ' Z podmíněné sdružené hustoty rozdělení (15) se vypočte podle definičního vzorce podmíněná ~třední hodnota složky '1]1

m~l = E{'1]l I '1]2 = A2:1b} (17)

Transformace (9) zajišťuje rovnost složek Sl a, '1]1' Proto jsou si rovny také podmíněné' střední hodnoty m. l a mw Velikost druhé složky m. vyplývá z trans­formace (9) a podmínky (2)

(18)

Složky mq a m.2 podmíněné střední hodnoty (7) jsou nejlepším odhadem chyb měření ve smyslu kri­téria (5), Předpoklad o rozdělení chyb měření je obec­nější ve srovnání s vyrovnávacím počtem, který před­pokládá rozdělení normální.

Y1Q

R2~ ____ -+ ______ ~~

Obr. 1. Postup určení opravy měr'ení s rovnoměrně rozloženými chybami

au~adzac. 260 ________________________________________________ ~3~1~(1~9~~)~č~.~10

PH klad 2

Nechť chyby měření z příkladu 1 mají rovnoměrné rozdělení li' hustotou p(xI ' x2 ) > O v oblasti ohraničené vobr. 1 obdélníkem RIR,fi,sRf, ' Vazební podmínku (3) znázorňuje přímka X2 = Xl + U:! -~. Transfor­mace (9) znamená přechod na souřad.iůcovou soustavu s osami YI a Y2' Osa Y2 je totožná. s osou X 2 a osa YI je natočena oproti Xl o 7tf4. Měřítko osy YI je znázorně­no polohou jednotky. V souřadnicích YI a Y2 má ná­hodná veličina 1]2 velikost U:! - ~ a náhodná veličina rJI má rovnoměrné rozdělení s nenulovou velikostí hustoty 1!t(YI I1]2 = U 2 -~) v intervalu YIP < YI < YIU'

Podmíněná střední hodnota veličiny '71 je Ledy dána aritmetickým průměrem souřadnic YIP a YIQ, které jsou shodné se souřadnicemi XIP a X1Q' Oprava údajú óbou čidel je dána souřadnicemi bodu S . Z grafického řešeni vyplývá, že oprava je závislá na. vzájemné poloze přímky představující vazební podmínku (3) a obdélníku RIR,fi,3R,.

4. Normálni rozděleni chyb měřeni

Sdružená. hustota rozdělení náhodných chyb měření s normálním rozdělením · je dána vzorcem

p(x) = [(27t)nIC lr~.exp(- (X-PE)T~-I(X-Pc»)

(19)

ve kterém p. označuje střední hodnotu a C varianční matici, ktel'á musí být symetrická. a pozitivně defi­nitní. Zkrácený zápis je e":"" Nn(pc C).

Pro zjednodušení výpočtu je vhodné použít sub­stituci

(20)

s horní trojúhelníkovou maticí H danou vztahem

C= HT.H (21)

K numerickému výpočtu rozkladu matice C lze použít Choleského algoritmus [3] . Rozdělení náhod­ného vektoru ~ je Nn (O, E). Pro novou proměnnou ~ má soustaV;l. rovnic (2) tvar

A (HT~ + f'.) = b (22)

Po substituci P = A.HT

q = b-A.pc (23)

má vazební podmínka formálně stejný tvar jako rov­nice (2)

P~= q (24)

Matici P je vhodné rozložit na bloky Pl a P2 stejným zpúsobem jako matici A v předcházející kapitole. Pro formulaci vztahú mezi bloky matice A a, matice P je třeba rozdělit do blokú ta,ké ma,tici H

H = [Hlll H12] O H22

(25)

Platí (26)

P 2 = A2Hf2 (27)

Podle obecného postupu výpočtu opravy z kapitoJy 3 je třeba, transformovat vektor ~ na vektOJ: 'I-

V transformační matici (10) jsou bloky matice A nahrazeny " bloky matice P. Rozdělení náhodného vektoru 'I je Nn(O, Q.QT). Vazební podmínka, pro vektor 'I je

(28)

" Po dosazení sdružených hustot rozdělení přísluš­ných normálnímu rozdělení za U a U2 do vzorcú (15) a (16) lze odvodit podmíněnou hustotu rozdělení Ul se střední hodnotou

(29)

Matice P. pT je čtvercová řádu m a má podle teorie matic stejnou hodnost jako P a A, a proto inverzní matice P. pT existuje. Podmíněné střední hodnoty složek ~l a '11 jsou si rovny ml;l = m1]l' SloŽka ml;2 je dána podobně jako (18) vzorcem

(30)

Po dosazení do (30) za vektor ml;l podle vzorce (29) plyne

(31)

Složky vektoru ml; dané vzorci (29) a (31) lze spojit

ml; = PT(P. PT)-l. q (32)

Vzhledem k substituci (23) je podmíněná střední hod­nota (7) a tím také oprava e* pozorovaných hodnot u·· dána vzorcem

s* = C.AT.(A.C.AT)-l.(b -A.p.) + f'. (33)

Numericky nejnáročnější operací při výpočtu (33) je inverze mat.ice řádu m. Rozklad (21), transformace (10) ll; (20) nejsou zapotřebí. Lze dokázat, že nezáleží na uspořádání sloupcú v matici A . .

Poznámka:

a) Vektor q ve vzorci (32) je násoben Moore-Penro­seovou inverzí

(34)

matice P, která minimalizuje euklidovskou délku vektoru ~ za podmínky minimální velikosti (P. ~ -_q)T. (P . ~ _ q) , která je podle (24) nulová. Výpočet oprav lze nahradit výpočtem vázaného extrému.

b) Pro p. = O a. diagonáliú tvar matice C je oprava podle (33) shodná s opravou vyrovnávacím počtem. Např. podle [5] představuje vektor (A.C.AT)-lb koreláty k, které se počítají z normálních rovnic

A.C.AT. k = b (35)

ve kterých je hlavní diagonála C sestavena. z převrá­cených hodnot vah popisujících přesnost měření slo­:lek vektoru e.

Příklad 3

Nechť chyby z příkladu 1 mají rozdělení e /'"OoJ N2 (p.u C), kde "

p. = [I-'Cl ], C = [l1.1 C12 ] (36)

1-'.2 C21 C22

Podle dříve odvozených vzorcú platí:

(37)

_3_1 --,(_19_8_8),---č_. _10 ________________________ __ 261 automatizace

(38)

Cl2 - C22 mE2 = . .q + flE2 Cu - 2,C12 + C22

(39)

Pro nezávislé chyby měření (C12 = O) bez systematické složky (PEl = flE2 = O) platí:

mEl = Cu

• (U1 - u2) (40) Cu + C22

mE2 = C22

,(U2- U l) (41) Cu + C22

Vzorce (40) a (41) lze odvodit také použitím postupu VÝp00tu známého z vyrovnávacího počtu.

5. Opravy parametrů sdružené hustoty rozděleni chyb měfenf

Při opakování měř'ení údajů čidel lze upřesňovat odhady parametrů sdružené hustoty rozdělení náhod­ných chyb. Vektor chyb e a vektor I? v soustavě lineárních rovnic (2) jsou v tomto případě náhodné funkce diskrétního času t = 1, 2, .. . Řada známých hodnot b (T) , T = 1, ... , t umožňuje vyhodnotit odhad střední hodnoty E{ b}, označený mb(t) , a odhad va­rianční matice var b , označený Cb(t). Užije-li se pro výpočet odhadů principu exponenciálního zapomínání (pravděpodobnost minulé hodnoty je IP-násobek pravděpodobnosti výskytů aktuální hodnoty, ° < IP, IP ~ 1), nabývají rekurzivní vzorce tvar

P(l) = ° P(t) = IPI(l - p(t - 1) + IP), t > 1

mb(t) = mb(t-l) .P(t) + b(t).(I-P(t»

Cb(t) = C b (t-l).P(t) + (b(t)-

- mb(t». (b (t) - mb(t))T. (1 - P(t» (42)

Mezní hodnotou rostoucího p(t) je parametl' IP· Vodítkem k volbě IP je například počet kroků, při kterém klesne pravděpodobnost výskytu příslušné zpožděné veličiny na polovinu:

IP I 0,999 0,995 0,99

i pro lPi = 0,5 . 692 138 68

0,98

34

0,95

14

0,9

7

Vektor b je podle (2) dán transformací vektoru e s transformační maticí A . Podle známých vět teorie pravděpodobnosti [4] platí

A . E{e} = E{b}

A . var e .AT = var b

(43)

(44)

Protože matice A je typu m X n s m < n, nelze složky vektoru flE a prvky matice C vypočítat řešením soustav rovnic (43) a (44). Jostliže PE a C budou považovány za náhodné veličiny s apriorně zvolenými sdruženými hustotami rozdělení P( fl E) a p(C), lze v časových oka­mžicích t opravovat flE a C na základě postupně se upl"esňujících vazebních podmínek

A. flE(t) = mb(t)

A . C(t) .AT = Cb(t)

(45)

(46)

postupem popsaným v kapitole 3. Vazební rovnice

(45) má tvar shodný s rovnicí (2) , a proto výpočtové vzorce se odvodí záměnou příslušných veličin. Rovni­ce (46) je soustavou s počtem m 2 rovnic s n 2 neznámými veličinami Clj, i, j = 1, ... , n. Tuto soustavu lze pře-psat do tvaru

(47)

kde tenzorový součin matic A 181 A je definován vzorcem

A 181 A = [ au A , ... , a1n A' ] , A = [a;!1 (48) am1 A , ... , amn A

a jednosloupcové matice [c]s a [Cb]s jsou utvořeny tak, že sloupce matice C = [Cl' ... , Cn], popř . Cb jsou uspořádány pod sebou do jednoho slOl,lpce

([C]s)T = (cr, ... , c;) (49)

Matice A ® A je typu m2 X n2 s hodností k( A ® A ) = = m2, tj. všechny rovnice soustavy (47), resp. (46) jsou navzájem nezávislé. Vzhledem k tomu, že matice C a Cb jsou symet.rické, lze omezit počet vazebních rovnic (47) na m(m + 1)/2 a počet neznámých veličin Clj , i = 1, ... , n , j = 1, ... , i, na n(n + 1)/2. Vyberou se rovnice (46) s pravými stranami Cbr3' r = 1, ... , m, 8 = 1, . . . , r, a položí se CJi = C;j'

6. Problém volby charakteristik normálně rozdělených chyb měřeni

Při jednorázové opravě měřených údajů je vhodné pro p(x ) ve tvaru (19) volit PE = O a prvky matice C s ohledem na známou třídu přesnosti jednotlivých typů čidel. Oprávněnost volby nulové střední hodnoty vyplývá z toho, že cejchování čidel není provedeno zá­měrně se systematickou chybou. Prvky na diagonále matice C jsou rozptyly chyb jednotlivých čidel. Na rozdíl od vyrovnávacího počtu je možné volit nenulo­vé prvky C mimo hlavní diagonálu jako důsledek závislosti některých čidel na společném zdroji poruch, kterým může být napříKlad teplota okolí čidel apod.

Přistoupí-li se při opakovaných měřeních k opra-. vám odhadu parametr:ů sdružené hustoty rozdělení chyb p( x ), je třeba zadat apriorní sdružené hustoty rozdělení P(flE) a p(C) , nejlépe ve tvaru příslušném normálnímu rozdělení. Pro p. ,....., Nn (O, c;.) je volena nulová střední hodnota. Varianční mat ici C~ lze sestavit podle dlouhodobé stálosti údajů čidel. Matici C~ je vhodné volit diagonální. Vypočteným opravám veličiny PE lze přisoudit význam systematických chyb čidel. Rozptyl údajů čidel určují prvky matice C . Uspořádá-li se k = n(n + 1)/2 teoreticky různých prvků matice C do sloupce [C]s, je rozdělení těchto prvků [C]. ,....., N k (pc, Ce) . Nenulové střední hodnoty pc je vhodné volit u prvků matice C z hlavní diago­nály. Varianční matici Ce je rovněž vhodné volit diagonální. Při zadání prvků této matice je třeba postupovat tak, aby při opravách rozptyly údajů čidel prakticky nedosáhly záporných hodnot. Roz­ptyly kovariancí z matice C musí být takové, aby tato matice zůstala po opravě pozitivně definitní. Pozitivní definitnost matice C je vhodné po opravě kontrolovat . Užitečnost opravy matice C je rovněž dána obtížností numerických výpočtů. Pro počet vazebních rovnic a zároveň řád invertovaných mat ic u opravy C platí:

eutomatlzace 262 __________________________ 3_'----'-('_9_88-.!-l_č_._'_o

m 1 2 3 4 5

m(m + 1)/2 1 3 6 10 15

Při klad 4

Tento příklad je pokračováním příkladů 1 a 3. Jsou v něm uvedeny výsledky dvou simulací oprav chyb měření a upřesňování parametrů jejich sdružené hustoty rozdělení. Bylo zvoleno:

8""'" N 2((0,5 , 2)T, diag (2, 1))

q; = 0,9

([C]s)T = (cll> C21 ' C22 )

[C]s ""'" N a((3, 0, l)T, diag (1 ; 0,09, 0,25))

1. simulace: f'. ,...." N 2 ((O, O)T, diag (2, 2)) 2. simulace: f'. ,...." N 2 ((0, O)T, diag (1, 3))

Po sérii 50 měření se simulovanými chybami byl vy­počten průměr druhých mocnin odchylek obou chyb měření od opravy vypočtené podle známých postupů z vyrovnávacího počtu {vzorce (40) a (41)) a podle postupu navržného v tomto článku (vzorce (38) a (39) včetně opravy parametrů hustoty rozdělení). vý­sledky jsou v tab. 1. Jak je patrno ze zadání, ve druhé Rimulaci byly předpokládány větší změny systematic­ké složky chyby Cl oproti chybě c2' Přesnost opravy využívající informace o základních statistických cha­rakteristikách změn vektoru b je větší než u klasické jednorázové opravy vypočtené vyrovnávacím počtem.

7. Závěr

V článku je odvozen obecný postup výpočtu opra­vy vzájemně závislých měřených údajů, který je konkretizován pro rovnoměrné a normální (vzorec (33)) rozdělení náhodných chyb měření. Výchozí

Tab. 1. Výsledky simulace oprav chyb měi'ení

Průměrná

I I Símulace

odchylka Vzorce

I 1. 2.

(40) 2,910 ( ' I - ,r)'

I (38) 2,171 1,360

(41) 2,910 ('2 - e~)2

I . (39) 1,876 1,256

předpoklady o statistických vlastnostech chyb měření jsou obecnější ve srovnání se známým výpočtem opravy závislého pozorování z vyrovnávacího počtu. Toho je využito k prúběžnému upřesňování odhadu párametrů sdružené hustoty rozdělení náhodných chyb měření na základě vyhodnocení opakovaných měření údajů čidel.

Výpočet oprav měřených údajúlze použít v automa­tizovaných systémech sběru dat. Absolutní velikost opravy, systematická složka a rozptyl chyb měření mohou sloužit ke kontrole stavu příslušného čidla.

I.iteratura

[1] HLAV.-\(~EK, Y. a kol.: Bilanční a simulační výpočty s loži­H'ch procesli chemické technologie. Praha, Aca.demia 1979.

[2] 'l'ÚMA, J.: Zpí'csnění údajú {:idel využitím r edundance v souboru informací o průběhu t.echnologického procesu. In: Sborník konfercncc Automatizace y hutnictví. Ostrava, Dům tcchniky ORVTS 1986, s. 26- 30.

[3 ] FADĚJEV, A. K. -- FADĚ.TEYOV.\, V. N.: výpočtové me­tody lineární a lgebry (rusky). Moskvu, Fizmatgiz 1960.

[4] ANDf<~L, J.: i\latcmatická "tatistika. Praha, SN'l'L jAJ-,FA 197 8.

[5] BAHTSCH, II. J. : Matematické vzorce. Praha, SNTL/ALFA 1971.

Došlo: 15. 9. 1987 L ektoroval: dor. Ing .. J . .Jenčík, CSc.