XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
239
Vybrané metody statistické regulace procesu pro autokorelovaná data NOSKIEVIČOVÁ, Darja Doc., Ing., CSc. Katedra kontroly a řízení jakosti, FMMI, Vysoká škola báňskáTechnická univerzita Ostrava, 17. Listopadu 15, 708 33 Ostrava-Poruba, Czech Republic
[email protected]
Abstrakt: V posledních letech se můžeme stále častěji setkat s kritikou klasických Shewhartových regulačních diagramů, ať už ze strany uživatelů či v odborných publikacích. Kritikové poukazují na to, že existuje velké množství procesů, kde je aplikace klasických metod nesprávná či nerealizovatelná. Praktikové na základě svých zkušeností často docházejí k závěru, že statistickou regulaci procesu (dále SPC) nelze v podmínkách jejich procesů aplikovat. Tento mylný závěr plyne z neznalosti předpokladů efektivní aplikace klasických Shewhartových regulačních diagramů, v lepším případě z absence ověřování splnění těchto předpokladů a z neznalosti netradičních metod SPC, které umožňují regulovat i procesy, u nichž nejsou uvedené předpoklady objektivně splněny. Dalším závažným problémem je nedostatečná znalost charakteru variability procesu. V první části příspěvku jsou definovány předpoklady efektivního praktického uplatnění klasických metod SPC a větší pozornost je věnována porušení předpokladu o vzájemné nezávislosti dat. Ve druhé části jsou pak diskutovány čtyři možnosti, jak problematiku autokorelace dat v rámci SPC řešit. Klíčová slova: autokorelovaná data, prodloužený kontrolní interval, statistika EWMA, ARIMA modely, dynamický diagram EWMA.
1 Předpoklady o datech a SPC Kromě obecného předpokladu efektivního využívání klasických metod SPC, tj. vyhovující způsobilosti měřicího systému, musí být splněna řada dalších předpokladů. Tyto další předpoklady lze rozdělit na dvě skupiny A) základní statistické předpoklady; B) ostatní. Mezi základní statistické předpoklady se řadí: 1. normální rozdělení znaku jakosti s konstantní střední hodnotou a konstantním rozptylem; 2. vzájemná statistická nezávislost hodnot znaku jakosti. Do ostatních předpokladů lze zařadit: 3. dostatečný počet dat; 4. citlivost na větší změny procesu; 5. sledování pouze 1 znaku jakosti na jednotce produktu.
2 Analýza předpokladu vzájemné nezávislosti dat Hlavním předpokladem efektivní aplikace klasických Shewhartových regulačních diagramů je vzájemná nezávislost hodnot sledovaného znaku jakosti. I velmi nízký stupeň vzájemné
240
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
závislosti hodnot sledovaného znaku jakosti (autokorelace) vyvolává selhání klasických Shewhartových regulačních diagramů. Selhání má podobu vysokého počtu zbytečných signálů. Tento jev není vůbec výjimečnou záležitostí v případě spojitých procesů, kde je autokorelace vyvolána velkou setrvačností procesu v čase. Stále častějším fenoménem se však vzájemná závislost dat stává i v podmínkách diskrétních procesů. Důvody lze spatřovat v automatizaci výrobních, zkušebních a kontrolních postupů, což umožňuje získat data z každého produktu (a je-li to potřebné, nejen jednoho znaku jakosti), tedy nejen z výběru n produktů odebraných z procesu po uplynutí určité doby Tv (kontrolního intervalu) od předchozího výběru, jak je obvyklé při realizaci sběru a záznamu dat při aplikaci klasických Shewhartových diagramů a jejich modifikací. (Pro kontrolní interval Tv nutno doplnit podmínku, že Tv » tv, kde tv je doba mezi 2 za sebou odebranými jednotkami ve výběru - viz obr. 1.) tv TV Obr. 1 Zobrazení principu tvorby logických podskupin při klasických regulačních diagramech Uvedený problém vzájemné se vyskytuje zejména u diskrétních procesů v krátkými výrobními cykly a s vysokou výrobní rychlostí, obecně řečeno u procesů s velmi krátkým intervalem mezi měřením a záznamem dvou po sobě jdoucích hodnot sledované veličiny. Existuje několik postupů, které vedou k odstranění autokorelace dat a umožňují využití SPC i v podmínkách vzájemné závislosti hodnot sledovaného znaku jakosti. Čtyři možnosti řešení uvedeného problému jsou diskutovány v další kapitole příspěvku.
3 SPC pro autokorelovaná data V následujících podkapitolách jsou předpokladu o nezávislosti dat:
podrobně rozebrány 4 varianty řešení
nesplnění
1.
metoda prodloužení kontrolního intervalu;
2.
postup s využitím aparátu modelování časových řad pomocí ARIMA modelování;
3.
aproximační postup založený na využití statistiky EWMA;
4.
dynamický EWMA diagram.
První tři varianty mají společný základ, který lze zobrazit následovně:
241
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
Vstup
i-tá část procesu
Výstup
Senzor
Vstup
(i+1)-tá část procesu
Výstup
Původní autokorelovaná data
x1, x2,...x t... FILTR
Nekorelovaná data xk,x2k,x3k…..nebo e1,e2,…..et…
Regulační diagram (klasický Shewhartův, EWMA, CUSUM)
xt, et
UCL CL LCL 1 2 3...
Číslo výběru (podskupiny)
k
Obr. 2 Zobrazení principu 1.-3. varianty odstranění autokorelace dat V případě 1. varianty představuje filtr (tj. prvek, který odstraní z dat korelační strukturu ) změna kontrolního schématu. Klasický Shewhartův regulační diagram pro individuální hodnoty je pak aplikován na každou k-tou hodnotu. U druhé varianty je filtrem vhodný stochastický model časové řady a regulační diagram se aplikuje na rezidua použitého modelu. Ve třetí variantě hraje úlohu filtru statistika EWMA. Vybraný regulační diagram se aplikuje na chyby této jednokrokové predikce.
3.1 Metoda prodloužení kontrolního intervalu Jak již bylo diskutováno, možnost získat údaje o sledovaném znaku či několika znacích jakosti z každé jednotky produktu přináší ze statistického hlediska překážku v podobě autokorelace dat. Nejjednodušším řešením je snížení frekvence statisticky zpracovávaných údajů. Bylo dokázáno, že s prodlužováním intervalu mezi 2 hodnotami sledované veličiny, které jsou statisticky zpracovávány, klesá míra korelace hodnot této sledované veličiny. V praxi to znamená, že zahrneme-li do zpracování dat hodnotu sledované veličiny z každé jednotky produktu, jsou tyto hodnoty za určitých podmínek vzájemně závislé. Zahrneme-li do zpracování dat každou k-tou hodnotu, pak s rostoucím k se snižuje míra korelace mezi daty, až do určité hodnoty k, kdy lze říci, že hodnoty sledované veličiny již spolu nekorelují. To lze postupně pro různé hodnoty k ověřovat pomocí grafu kolerace nebo pomocí autokorelační funkce. Uvedené řešení je velmi jednoduché, avšak jeho aplikace znamená, že není využita celá informace, získaná o procesu ze všech naměřených hodnot. Jestliže např. zpracujeme každou 10. hodnotu, nevyužijeme 90% informace, která je k dispozici. Další varianty řešení autokorelace představují složitější, ale vhodnější řešení.
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
242
3.2 SPC s využitím ARIMA modelů Postatou tohoto řešení je nalezení vhodného modelu časové řady a aplikace regulačního diagramu na rezidua modelu (odchylky skutečně naměřené hodnoty od hodnoty vypočtené dle modelu). Jako nejvhodnější se jeví Box-Jenkinsovy stochastické ARIMA modely (Autoregressive Integrated Moving Average) [8]. Box-Jenkinsonova metodika představuje moderní koncepci analýzy stacionárních a nestacionárních časových řad, založenou na teorii pravděpodobnosti. Obecný tvar modelu ARIMA(p,d,q) je následující:
Φp(B)∇ d xt = Θq(Β) εt
(1)
kde Φp(B) = (1- φ1Β - φ2Β2 - ... - φpΒp) je autoregresní polynom p-tého řádu, Θq(Β) = (1 - θ1Β - θ2Β2 - ... - θqΒq) je polynom klouzavých průměrů q-tého řádu, ∇ je operátor zpětné diference (tento prvek se zavádí v případě, že modelovaný proces vykazuje nestacionaritu), d je řád diference, B je operátor zpětného posunu ( B.xt = xt-1), φ1, φ2,..., φp jsou parametry autoregresního modelu, θ1, θ2,...θq jsou parametry modelu klouzavých průměrů. εt je proměnná, které se říká bílý šum a představuje nepredikovatelnou fluktuaci v datech. Má normální rozdělení se střední hodnotou rovnou nule a konstantním rozptylem σp2 a její hodnoty jsou nekorelované. Je-li ˆx t odhad empirické hodnoty xt získaný pomocí vhodně zvoleného ARIMA modelu, pak rezidua tohoto modelu et = xt - xˆ t se budou chovat jako nezávislé náhodné proměnné pocházející z normálního rozdělení. Nejčastěji se v praxi používají následující ARIMA modely. Uvažujme model xt = ξ + φxt −1 + ε t
(2)
kde ξ a ϕ(-1< φ<1) jsou neznámé konstanty a εt je normálně rozdělená a nezávislá veličina se střední hodnotou rovnou nule a směrodatnou odchylkou σ. Tento model se nazývá autoregresní model 1. řádu a označuje se AR(1). Hodnoty sledovaného znaku jakosti, které jsou navzájem posunuté o k časových period (xt a xt-k), mají korelační koeficient φ k . To znamená, že autokorelační funkce ACF by měla exponenciálně klesat. Rozšíříme-li rovnici (2) do tvaru xt = ξ + φ 1 xt −1 + φ 2 xt − 2 + ε t ,
(3)
dostáváme rovnici autoregresního modelu druhého řádu AR(2). Obecně je v autoregresních modelech AR(p) proměnná xt přímo závislá na předchozích hodnotách xt-1, xt-2, atd. Jestliže modelujeme závislost dat pomocí náhodné složky ε t , pak dostáváme modely klouzavých průměrů MA(q). Model klouzavých průměrů 1. řádu má rovnici: xt = µ + ε t − θε t −1
(4)
V tomto modelu je nenulová korelace pouze mezi dvěma po sobě jdoucími hodnotami xt a xt-1 a lze ji vyjádřit následovně: ρ 1 = −θ /( 1 + θ 2 ) . Tomu odpovídá tvar autokorelační funkce ACF (Arlt, 1999).
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
243
Pro modelování praktických úloh se často hodí složený model obsahující jak autoregresní složku, tak složku klouzavých průměrů. Tento model se obecně označuje ARMA(p,q). Model ARMA 1. řádu, tj. ARMA(1,1) má rovnici: xt = ξ + φxt −1 + ε t - θε t -1
(5)
Tento model je často vhodný pro chemické a jiné spojité procesy, kde modelem AR(1) lze velmi dobře modelovat mnohé znaky jakosti. Náhodnou složkou modelu jsou pak popsány chyby měření, o kterých předpokládáme, že jsou náhodné a nekorelované. V ARMA modelech se předpokládá stacionarita procesu, tzn., že hodnoty sledovaného znaku jakosti se pohybují kolem stabilní střední hodnoty. Avšak často se v praxi objevují procesy ( např. v chemickém průmyslu, kde sledovaný znak jakosti xt je výstupní veličinou, který není žádnou regulací udržován na cílové hodnotě), kde hodnoty sledované veličiny „utíkají“. Pak je vhodné modelovat procesy pomocí vhodného modelu s operátorem zpětné diference ∇ , např. modelem ARIMA (0,1,1), jehož rovnice je: (6) xt = xt −1 + ε t − θε t −1 . Modely ARIMA vyhlížejí odlišně od Shewhartova modelu ( xt = µ + ε t pro t=1,2…). Jestliže však do rovnice (2) dosadíme za ϕ = 0 nebo v rovnici (4) za θ = 0 , dostaneme Shewhartův model procesu. Dalším důležitým krokem při využití ARIMA modelů pro SPC je volba vhodného regulačního diagramu. Je-li testováním reziduí zvoleného ARIMA modelu prokázáno, že jsou nekorelovaná a pocházejí z normálního rozdělení, je možné ověřit pomocí reziduí, zda proces je či není statisticky zvládnutý (zda působí nebo nepůsobí vymezitelné příčiny). Protože rozsah výběru n = 1 (původní empirické hodnoty xt byly zjištěny u každé vyráběné jednotky), nabízí se na prvním místě dvojice regulačních diagramů pro individuální hodnoty a klouzavé rozpětí. Střední přímka CL a horní a dolní regulační meze UCL a LCL se u diagramu pro individuální hodnoty stanoví ze vztahů: CL = e (≅ 0),
kde
(7)
UCL = e +
3 R kl , 1.128
(8)
LCL = e −
3 Rkl , 1.128
(9)
e je průměrná hodnota reziduí, Rkl je průměrné klouzavé rozpětí .
Střední přímka CL a regulační meze UCL a LCL v diagramu pro klouzavé rozpětí se určí takto: CL = Rkl ,
(10)
UCL = 3.267 Rkl ,
(11)
LCL = 0. (12) Chceme-li zvýšit citlivost regulačních diagramů reziduí z ARIMA modelů na menší odchylky, doporučuje se použít oboustranný regulační diagram CUSUM (Cumulative Sums) s rozhodovacím intervalem ±H nebo klasický diagram EWMA, oba aplikované na rezidua (Montgomery, Friedman, 1989). Návrh oboustranného diagramu CUSUM spočívá ve stanovení rozhodovacího intervalu ±H a parametru K . Do oboustranného regulačního diagramu CUSUM s rozhodovacími mezemi
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
244
+H a –H, aplikovaného na rezidua, se zaznamenávají hodnoty kumulovaných součtů vypočtených dle vzorců: S+t = max [0, S+t-1 + (et - e - K)],
(13)
(14) S-t = min [0, S-t-1 + (et - e + K)]. Alternativním diagramem k výše uvedenému diagramu CUSUM je EWMA diagram pro rezidua. CL, UCL a LCL v tomto diagramu vypočteme dle vztahů:
UCL = e + ∆
LCL = e - ∆
CL = e (≅ 0),
(15)
R kl
λ 2t , 1 − ( 1 − λ ) ( 2 − λ )
(16)
λ 2t . 1 − ( 1 − λ ) (2− λ)
(17)
1.128 R kl 1.128
λ (0 < λ < 1) je parametr zapomínání a ∆ je vzdálenost regulačních mezí od střední přímky v počtu směrodatných odchylek. Obvykle se doporučuje volit malé hodnoty λ (např. (0.05 ≤ λ ≤ 0.2) a buď ∆ = 2.5 nebo 3 (Montgomery, 2001). Sleduje-li se na jednom produktu m znaků jakosti najednou, je možné aplikovat na rezidua z m ARIMA modelů Hotellingův T2 diagram nebo CUSUM či EWMA diagram pro vícerozměrné proměnné (Montgomery, Friedman, 1989). Ukázky regulačního diagramu pro individuální hodnoty a klasického EWMA diagramu, aplikovaných na rezidua ARIMA modelu, jsou na obr. 3a a 3b.
E
Obr. 3a
Regulační diagram pro individuální hodnoty aplikovaný na rezidua zvoleného ARIMA modelu
Obr. 3b EWMA regulační diagram pro rezidua zvoleného ARIMA modelu
3.3 Aproximační postup založený na využití statistiky EWMA Modelování procesu pomocí ARIMA modelů není z praktického hlediska příliš jednoduchou záležitostí, i když dnes každý kvalitní statistický softwarový balík ARIMA modely obsahuje. Montgomery a Mastrangelo (1991) vytvořili návrh postupu, který je aproximací přesnějšího ARIMA modelování. Tento aproximační postup je založen na použití statistiky EWMA. V postupu je využito faktu, že EWMA může být za určitých podmínek použito i pro
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
245
autokorelovaná data. Předpokládejme, že proces můžeme modelovat pomocí modelu ARIMA(0,1,1), popsaného rovnicí xt = xt −1 + ε t − θε t −1 .
(18)
EWMA s parametrem λ = 1 − θ je optimální jednokrokovou predikcí pro modelovaný proces. Konkrétně, jestliže ˆxt +1 je předpovědí hodnoty sledované veličiny pro časový okamžik t+1, zrealizovanou v časovém okamžiku t, pak ˆxt +1 = EWMAt = λ .xt + ( 1 − λ )EWMAt −1 .
(19)
Chyba této predikce v časovém okamžiku t se stanoví ze vztahu: et = xt − ˆxt ,
(20)
kde xt ……….je hodnota sledované veličiny v časovém okamžiku t, ˆxt ………je odhad hodnoty sledované veličiny v časovém okamžiku t provedený v časovém okamžiku t-1. Hodnoty chyby predikce mají normální rozdělení s nulovou střední hodnotou a nejsou korelovány. Statistická regulace procesu je pak realizována tak, že se regulační diagram pro individuální hodnoty aplikuje na chyby jednokrokové predikce et (tj . na rezidua modelu EWMA). Parametr λ by měl být stanoven metodou nejmenších čtverců chyb et. Uvedený postup lze použít i na procesy, pro které by byl vhodnější jiný model než diskutovaný model ARIMA(0,1,1). Obecně platí, že jestliže hodnoty sledovaného znaku jakosti jsou pozitivně korelovány a střední hodnota procesu se mění pomalu („Slow Drift“), pak EWMA s vhodnou hodnotou parametru λ poskytuje výbornou jednokrokovou predikci. (Montgomery, 2001). Společně s regulačním diagramem reziduí by měl být veden graf jednotlivých skutečně naměřených hodnot, kam by se souběžně měly zaznamenávat hodnoty EWMA (Montgomery, 2001).Tímto způsobem je informace o statistické zvládnutosti procesu obsažená v diagramu reziduí doplněna o vizualizaci dynamiky procesu a přesnosti odhadů.
3.4 Dynamický EWMA diagram V tomto případě obsahuje informaci o statistické zvládnutosti procesu i o dynamice procesu jediný diagram (na rozdíl od předchozí varianty, kde je doporučeno vést jak regulační diagram pro rezidua EWMA modelu, tak grafy skutečně naměřených hodnot sledované veličiny a odhadů hodnot této veličiny - Montgomery & Mastrangelo, 1991) Uvedený dynamický diagram EWMA je vhodný pro jmenované procesy tehdy, když hodnoty sledované veličiny vykazují pozitivní autokorelaci a proces má nekonstantní střední hodnotu s pomalými změnami. Překročení regulačních mezí způsobí v tomto diagramu pouze náhlá změna střední hodnoty, malé změny procesu diagram toleruje. Dynamický diagram EWMA poskytuje tedy jak informaci o statistické zvládnutosti procesu, tak o jeho dynamice. Princip dynamického EWMA diagramu lze popsat následovně: Je-li statistika EWMA vhodnou jednokrokovou predikcí, pak lze hodnoty EWMAt použít jako hodnoty střední čáry CL v časovém okamžiku t+1. Horní a dolní regulační mez lze pak stanovit ze vztahů UCLt +1 = EWMAt + 3σ p ,
(21)
LCLt +1 = EWMAt − 3σ p .
(22)
XXVIII. ASR '2003 Seminar, Instruments and Control, Ostrava, May 6, 2003
246
S těmito regulačními mezemi jsou potom pro posouzení statistické zvládnutosti procesu porovnávány naměřené hodnoty sledovaného znaku jakosti xt+1. V dynamickém diagramu EWMA jsou jak střední čára (nejde o přímku jako u klasických diagramů), tak obě regulační meze dynamické. Určíme je dle vztahů CLt +1 = ˆy t +1 = EWMAt ,
(23)
LCLt +1 = ˆy t +1 − u 1−α / 2 .σˆ p = EWMAt − u 1−α / 2 .σˆ p ,
(24)
UCLt +1 = ˆy t +1 + u1−α / 2 .σˆ p = EWMAt + u1−α / 2 .σˆ p ,
(25)
kde u 1−α / 2 je kvantil normovaného normálního rozdělení.
UCL
xk
CL LCL
č.naměřené hodnoty
Obr. 4 Dynamický diagram EWMA
4 Závěr Setkáme-li se při praktické aplikaci SPC v problémem autokorelovanosti dat, můžeme využít jedné ze čtyř metod diskutovaných v tomto příspěvku. Příspěvek obsahuje jak popis jednotlivých metod, tak jejich stručné srovnání.
5 Literatura ARLT, J.1999. Moderní metody modelování ekonomických časových řad. Praha, Grada Publishing, 1999. ISBN 80-7169-539-4. MONTGOMERY D.,C. & FRIEDMAN,D.J.1989. Statistical Process Control in a ComputerIntegrated Manufacturing Environment. In: Statistical process control in automated manufacturing. New York, Marcel Dekker, Inc., 1989, pp. 67-89. MONTGOMERY, D.C. & MASTRANGELO, CH. M.1991. Some Statistical Control Methods for Autocorrelated Data. Journal of Quality Technology, 1991, sv. 23, č. 3, s. 179-193. MONTGOMERY, D.C.: Introduciton to Statistical Quality Control.2001. J.Wiley & Sons, New York, 2001. 796 s. ISBN 0-471-31648-2.