A sztochasztikus approximáció konvergenciájáról
 
 
2019 június
 
 
 
 
ABSTRACT
A sztochasztikus approximáció egy paraméterbecslő, a sztochasztikus folyamatok modellezésénél gyakran alkalmazott rekurzív algoritmus, melynek konvergencia sebessége 1/n szerű, n = 1,2,3,... . Az algoritmus sebességét segéd algoritmusokkal gyorsítjuk, a konvergencia gyorsításának következménye sok esetben torzított becslés. Relatíve kis torzítás általában elfogadható, pl. a modell fokszámának csökkentésekor, mert gyorsító hatású. A dolgozatban a predikciós hibafolyamat elemzésén alapuló algoritmust, és a hibák előjelein alapuló konvergencia sebességet szabályozó algoritmust mutatunk be, amely a predikciós hibafolyamat zérus várható értékét és előjel szimmetriáját írja elő. / On the convergence of stochastic approximation: The stochastic approximation is a recursive algorithm with a convergence rate of 1/n, n = 1,2,3,... , its speed is usually accelerated by an auxiliary algorithm. The consequence of speeding up convergence is often a biased estimate, but relatively small biases are considered acceptable for applications, even at the cost of reducing the number of degrees in the model, because of the speed-up effect. In this paper, we present a convergence rate control algorithm based on the analysis of the prediction error processes, the sign of the errors, which imposes a zero expected value and sign symmetry on the prediction error process. 
 
BEVEZETÉS
 
A sztochasztikus approximáció algoritmusa (STA) egy heurisztikus rekurzív paraméter becslő algoritmus, a matematikai modellezésnél gyakran alkalmazzuk. Az STA  a  relatív gyakoriság rekurzív becslésének rekurzív alakjából is származható: ha n = 1,2,3,4,...
 
 
Aátl.n = 1/n ΣAi = 1/n {(n-1) Aátl.n-1 + An } = (1- 1/n) Aátl.n-1 + A/ n = Aátl.n-1 + 1/n (An - Aátl.n-1)
 
Azaz egy Θn  paramétervektor értékének új becslését úgy számítjuk, hogy az érték előző lépésben számított becsléséhez hozzáadjuk az új  lépésben kapott hibájának 1/n -szeresét. A hibafolyamat matematikailag 0 várható értékű, véges szórású folyamat. Ha az STA hibafolyamat az output mért és az előző lépésben becsült értékeinek különbségeként határozzuk meg, akkor predikciós hibafolyamat, amely aszimptotikusan független az yn becslési folyamattól, általában jól teljesülő feltételek esetén, ld. a Kálmán-szürőt. (https://hu.wikipedia.org/wiki/K%C3%A1lm%C3%A1n-sz%C5%B1r%C5%91). A predikciós hibafolyamat és az yn -el jelölt becslési folyamat  korrelálatlanságának (aszimptotikus függetlensége) esetén a paramétervektor becslése legkisebb négyzetes értelemben optimális. Konvergenciája 1/n → 1/(n+ const.) helyettesítés esetén bizonyított. Megj.: exponenciális felejtés esetén (1 - 1/ν) Aátl.n-1 + 1/ν (An - Aátl.n-1) alakban számítjuk, ahol ν szokásos értéke 0.99 - 0.999. Az exponenciális felejtés jó kezdőértékek aktualizálására alkalmas algoritmus. Jó a kezdőérték, ha predikciós hibák becsült számtani közepe nulla.
 
A tekintett STA algoritmusoknál a hiba predikciós hibafolyamat, összetevői: becslési zaj, mérési zaj és rendszer zaj, a becslési haj négyzetösszegét minimalizáljuk.  A becsülendő paraméter vektort, Θ-t,  rekurzív alakban úgy módosítjuk,  hogy kiszámítjuk a módosítás irányát egy gn-el jelölt gradiens vektor alapján. A ggradiens vektor a predikciós hibák négyzet-összegének** a keresett paramétervektor szerint számított  gradiense (https://web.archive.org/web/20151112004321/http://www.sztaki.hu/~szcsaba/talks/lecture1.pdf). 
 
Az STA egy paramétervektor becslő algoritmus: a rekurzív relatív gyakoriság és a gradiensen alapuló (szélsőérték kereső) algoritmusok heurisztikus kombinációja. (A kimenőjelre vonatkozóan fel szoktuk tenni, hogy nulla várható értékű: ha nem ismert a várható érték akkor pl. az  yn,mért yn-1,mért  stacionáriusnak gondolt folyamatra számítjuk a prediktort**.  Az yn,mért folyamat lehet diszkrét is, itt folytonosnak tekintjük. A gradiens vektornak nNemlineáris komponensei is lehetnek, ami egy általánosítási lehetőség, továbbá a vektor output esete is egy lehetséges általánosítás. Nevezik az STA-t parametrikus tanuló algoritmusnak is (https://en.wikipedia.org/wiki/Naive_Bayes_classifier).
 
A gradiens alapú keresés (https://en.wikipedia.org/wiki/Stochastic_gradient_descent) esetén, ha a becsült yn (kimenet, output jelentésű) függvény a gradiensvektor lineáris függvénye, akkor  ynΘn T gn , ahol T transzponálást jelöl.  A gn gradiens vektorkomponensei  lineáris esetben a kimenet korábbi becslései, a hibafolyamat korábbi értékei, segédváltozó (az input standardizált) értékei, de lehetnek nemlineáris komponensek is:  a becslés Θ paraméterben lineáris. A becsült paramétervektor, egyben a gradiens vektor dimenziójának növelésével a becslés konvergenciája, azaz a tanulási arány gyorsan romlik, pl a hibafolyamat oszcillálva divergálhat is. Segédváltozó nélküli esetben egy fehér zaj input folyamatot és egy racionális törtfüggvénnyel, transzferfüggvénnyel leírható rendszert tételezünk fel. Az fehérhaj folyamat szórásnégyzeténak nagysága a folyamatra vonatkózó ismereteink bizonytalanságával arányos.  A  transzferfüggvény ismeretlen paramétereivel írjuk fel a gradiens vektort. 
 
A gradiens szerinti, "(Cauchy) a legnagyobb függvényérték-csökkenés iránya szerinti minimum keresési módszer -tehát  keresési módszer az STA-tól független - a legmeredekebb esés módszere. Iterációnként az aktuális pontban kiszámoljuk a deriváltat, ami a gn gradiens, majd vele ellentétes irányban, egy lépésben vagy kereséssel meghatározzuk az ebben az irányban elérhető LS, azaz legkisebb négyzetes célfüggvény minimumot. A minimumpont a következő iteráció kezdőpontja. " Az eredeti STA szerint a paramétervektor módosítása minden lépésben: a tanulási arány, a hibafolyamat és a gradiens szorzata (elvi alak skalár output yn,mért esetén):
 
                                                                                   ΘΘn-1 + 1/n  (y- yn,mért) gn
 
 
A sztochasztikus approximáció egy kb. 70 éves algoritmus család, sok változata ismert (https://en.wikipedia.org/wiki/Stochastic_approximation). Itt a gradiens egy normalizált változata fog szerepelni, ezért nem divergál a nagy gradiensű esetekben: 
 
                                                                        ΘΘn-1 + 1/n  (y- yn,mért(Σ gi2)-1/2 g,  
 
feltéve, hogy gi komponensei nem egyenlőek nullával (könnyen átírható mátrix alakba). 
Az algoritmus a gradiens nagy értékeire -  a (Σ gi2)-1/2 , és i = 1,2,..., dim g = dim   Θn) tényezők következtében- kevésbé érzékeny mint az eredeti STA algoritmus. A kis gradiensű eset numerikusan egyszerűen kezelhető. Nagy gradiensű esetekben célszerű lehet a gradienst a mozgó átlagával helyettesíteni. 
A továbbiakban a tanulási arányt és a hibafolyamatot vizsgáljuk. A hibafolyamat egyik összetevőjének oka: Θn  eltérése a valódi értékétől, azaz a ΔΘnT gn.  hibafolyamat. A hibafolyamat az ismeretlen ΔΘ -el lineárisan arányos. Ha 1/n gyorsabban tart nullához, mint ΔΘ, az torzított becslést okoz, ami a hibafolyamat tulajadonságaiból számított cn skalár szorzóval próbálunk elkerülni. Az 1/n skalár szorzó neve tanulási sebesség vagy arány (az angol terminológia szerint learning rate),  O (1/n) szerint csökken, ahol O ordót jelöl. Az állítás visszacsatolt lineáris irányítási (closed-loop) rendszerekre is érvényes feltéve, hogy egy lépés késleltetés van a hurokban. A bizonyítás a becslési hibára felírt Steiner-formula alapján történhet.(https://en.wikipedia.org/wiki/Bias_of_an_estimator).
  
A tanulási sebesség és a hibafolyamat konvergenciája
A konvergencia sebességét lehet javítani O (1/n) -hez viszonyítva, de ennek torzított becslés lehet a következménye. Szokásosnak nevezhető kompromisszum, hogy kevesebb paramétert becsülve - azaz a modell fokszámát csökkentve- gyorsul a konvergencia. Amennyiben 1/n -nél lényegesen gyorsabb a csökkenés, az le is állíthatja a rekurzív becslést, torzított becslést eredményez. Pl. ha jó kezdő érték meghatározásához minden lépésben külön iterációval nulláznánk a hibát, akkor zajkövető beccslést kapunk.(https://en.wikipedia.org/wiki/Bias%E2%80%93variance_tradeoff) A másik határesetben, a  paraméter komponensek abszolut értékeihez viszonyítva  nagy, de korlátos hibák esetén: a becslések oszcillálnak, majd divergálnak, ami a lépésnagyság mesterséges, legalább 1/n szerű csökkenésével megelőzhető. (https://pp.bme.hu/ee/article/view/4974, Bencsik, I. (1974) “CONVERGENCE RATE ACCELERATION OF SUCCESSIVE APPROXIMATION ALGORITHMS IN REAL-TIME CASE ”, Periodica Polytechnica Electrical Engineering (Archives), 18(1), pp. 99-103.) Gyakori eset, hogy kis torzított becslést elfogadva gyorsul a konvergencia.
 
Csak szuboptimális STA algoritmusok léteznek, de a szuboptimális tanulást is optimális rész stratégiák alkotják. A hiányos a priori ismereteink miatt a szuboptimális STA-nál: hibafolyamat 1/n szerű csökkenésével és nulla várható értékének biztosításával lehet elérni a konvergenciát. 
Továbbiakban csak a hibafolyamatot vizsgáljuk. A hibafolyamat minden olyan folyamatot tartalmaz, amit nem építettünk be a modellbe: ynΘnT g-be. Ezek a rendszerzaj, a mérési zaj és a "becslési zaj" (ez utóbbi ΔΘn -el kapcsolatos, tartalmazza a modellezéssel kapcsolatos pontatlanságokat, pl. a modell relative alacsony fokszáma következtében***). A paraméterbecslési módszereket meg lehet különböztetni a zajra vonatkozó feltételezések szerint is. Az STA esetén véges szórásnégyzetű, független, nulla várható értékű  fehér zajokat, és predikciós típusú hiba folyamatot feltételeztünk.
Gyakran alkalmazott módszer - és éppen az STA algoritmusokat és az y  = ΘnT g alakú modelleket használják a célra- egy korábbi becslés hibafolyamatának fehérítése, azaz a hibafolyamat korrelált részének becslése és levonása (angol terminológia: whitening the residual process). A Θn  - Θ változó zérus várható értéke biztosítandó (https://en.wikipedia.org/wiki/Fisher_information).
 
Egy, a mérési zajjal kapcsolatos probléma: egy jel valódi értékének becslése ismert szórású független mérési zaj esetén. Pl. ha a jel mérési zajos fehér zaj lenne és a szórásnégyzete adott lenne, akkor a valódi érték becslése:  a fehér zaj mért pillanatnyi értéke szorzandó a
 
(fehér zaj szórásnégyzete - a mérési zaj szórásnégyzete)/fehér zaj szórásnégyzete
 
hányadossal, a Kálmán-szűrő alapján. A rendszer zajjal - kisebb részben a mérési zajjal- kapcsolatos probléma a becslések zaj-érzékenysége: a jó kezdőértékek  becslésénél egy tanuló minta segíthet kezelni a problémát.  (A Kálmán-szűrő alapján belátható, hogy a kimenet és a paraméterek legkisebb négyzetes értelemben valóban optimális becsléshez, az u.n. extended Kalman-filterhez, a fokszám és a rendszerzaj, azaz a fehér zaj szórásnégyzetének és a mérési zaj szórásnégyzetének ismerte is szükséges. A becslés egybe esik a Bayes becsléssel, és a MLH becslésnél a filtráció járulékos beépítése szükséges.)
 
A hibafolyamat zérus várható értékét , a szuboptimális konvergencia sebességet a lépés nagyságának speciális megválasztásával biztosítjuk: egy cn -el jelölt, minden lépésben újra számolt skalár szorzó alkalmas megválasztásával: nem csak azt követeljük meg, hogy a hibafolyamat 1/n szerűen csökkenjen, hanem azt is hogy a hibafolyamat előjele a lépések ≈50% -ában különbözzön.( Előjelpróba: Vincze István: Matematikai Statisztika, Műszaki Könyvkiadó, 1968, 160. old.) A tanulási arányt ekkor cn/n  jellemzi. A bemenő és általában a jelek változékonysága (angol terminológia persistenly exciting) azért fontos, mert a jelek kellő változékonysága nélkül nem módosul becslés a szükséges mértékben. Ezt a változékonyságot részben biztosíthatja, hogy a lépésnagyság cn skalár szorzóját minden lépésben úgy számítjuk*, hogy a hiba - egy egységnyi skalár szorzóval tett próbalépés után- előjelet váltson és egyidejűleg 1/n szerűen csökkenjen abszolút értékben: ez a szuboptimális tanulási sebesség STA algoritmusok esetén. A legmeredekebb esés módszerét tehát a zérus várható értákű hibafolyamat előjelére is alkalmazzuk.*  A cn › 0 skalárral 1/n -él gyorsabb konvergencia is beállítható, de a gyorsabb tanulás ára a torzított becslés. 
KÖVETKEZTETÉSEK
A gyorsító algoritmusok - nagy számban léteznek- okozta torzítás csökkentésének módját vizsgáltuk.  A  legmeredekebb esés módszerét használó algoritmusok szuboptimális algoritmusok abban az értelemben, hogy a további gyorsítás már növekvő torzítást okoz. A sztochasztikus approximáció esetén a gyorsító módszer az 1/n szerűen csökkenő hibafolyamat előjel szimmetriáját írja elő. 
 
 
A ZÉRUS VÁRHATÓ ÉRTÉKŰ PREDIKCIÓS HIBA FOLYAMATRÓL 
 
A becslési, a rendszer és a mérési zaj folyamatokról feltesszük, hogy nulla várható értékű, fehér Gauss-folyamatok és korrelálatlanok a vizsgált és modellezendő folyamattal. A feltételezett zaj modellek befolyásolják a modellezést. A predikciós hibafolyamat kiértékelésére egy rövid, az utolsó 20 elemet tartalmazó mintát fogunk használni. A hibák előjeleinek szimmetriáját tesztelve található egy olyan  minden lépésben újra számított c› 0  skalár, amellyel biztosítjuk, hogy a hiba abszolut értékben O(1/n) szerint csökkenjen és a hibák előjelei azonos számban forduljanak elő. Ez utóbbit egy indikátor változó segítségével ellenőrizzük: meghatározzuk az azonos előjelek számának várható értékét. 
 
Állítás: Valamely b alapú számrendszerben, ha 1/b a számjegyek előfordulásának valószínűsége, jelölje P(X=k) egy k hosszúságú független kísérletekkel kapott tetszőleges mintázat valószínűségét megszámlálhatóan végtelen sokaság esetén. ( Pl. a pi szám bináris alakjának megszámlálható része. A mintázat fogalma a jegyek különbözőségével kapcsolatos.) Ekkor általános b-re:  
 
P(X=k) = p (1-p)k-1 (b-1)/b 1/b k-1 = (b-1)/bk  , for k = 1, 2, 3, ...,∞ 
 
ahol (b-1)/b a geometriai eloszlás paramétere, egyben a szinglik valószínűsége.  Az X valószínűségi változó várható értéke b/(b-1), ami a b-1)/b paraméter reciproka, amikor az eloszlás entrópiája maximum az exponenciális eloszlás családon belül. Az X változó szórása b/(b-1)2 -al egyenlő. b=2 esetben a várható érték 2A szinglik nélküli várható érték  (2b-1)/(b-1) = 3, , ekkor k = 2,3,4,.. , -al egyenlő.
Véges, kmax  nagyságú minta esetén a (b -1) Σk b-k valószínűségek összege:
(b -1) Σk b-k  = 1- b-kmax   ≈   1/2,  k = 2,3.., amit max 20 -nak választottunk.
 
Tehát a szinglik nélkül számított P(X=k) valószínűségek összege egy 20 elemet tartalmazó ablakban jó közelítéssel 1/2, és a szinglik előfordulásának valószínűsége is 1/2.
Az állítás alapján a 20 elemet tartalmazó ablakban a hibák előjeleinek előfordulási gyakoriságai értékelhetőek, és cn értékei számíthatóak az STA algoritmusban úgy, hogy az előjelek szimmetrikusan forduljanak elő. Egy további lehetséges módja lehet a lépésnagyságot módosító cn számításnak, hogy az yn,measuredy hibafolyamat abszolut értékére előírjuk, hogy az 1/n és 1/n1/2 közötti intervallumban csökkenjen, és minden lépésben előjelet váltson.

**

Elméleti kritérium is található az yn,measuredy hibafolyamat négyzetes várható értékre vonatkozó M(y2) = M2(y) + D2 (y) összefüggés alapján, feltéve, hogy M(y2) létezik:

                M2(y) /  M(y2) + D2 (y) / M(y2) = 1, 

tehát  D2 (y) / M(y2)  ≤ 1.  A  D2 (y) / M(y2)   relatív szórásnégyzet, értéke jól becsülhető, alkalmas legkisebb négyzetes mérőszám a becslések összehasonlítására. (https://en.wikipedia.org/wiki/Bias_of_an_estimator).

A  D2 (y) / M(y2) arány, a relatív szórásnégyzet hasznos kritérium az azonos és nem nulla várható értékű becslések kiértékelésénél, mert a Csebisev egyenlőtlenségben (https://en.wikipedia.org/wiki/Chebyshev%27s_inequality) ε2 helyére ε2 M (y2) -t helyettesítve, az egyenlőtlenség jobb oldala D2 (y) / ε2 M(y2) alakú lesz:

  P y - M (y)|  ≤ ε M1/2  (y2 } ≤   ε-2 (1 -  D2 (y) /  M (y2))            .

(Ld.: https://bencsik.rs3.hu/component/content/category/173-valoszinusegeloszlasok-informacio-tartalmanak-ertekelese-a-mernoeki-gyakorlatban.html?layout=blog&Itemid=101)

***

A modell fokszámának növelése "költséggel" jár, egy használható megoldás: jelölje az N a dim Θ -t, a paramétervektor és a modell dimenziószámát és D2N (y) / M(y2) a kritérium értékét. Ekkor érdemes agy N + 1 dimenziós modellt számítani, ha

                                                                                                 (N + 1) D2N+1 (y) < N D2N (y).

Tehát a magasabb dimenziós modell átlagos hibája még a dimenziójával szorozva is kisebb...