Elektrotehniški vestnik 76(4): 240-245, 2009 Electrotechnical Review, Ljubljana, Slovenija Modeliranje časovnih vrst z metodami teorije informacij Marko Bratina1, Andrej Dobnikar2, Uroš Lotrič2 1 Savatech, d.o.o., Škofjeloška 6, Kranj, Slovenija 2 Univerza v Ljubljani, Fakulteta za računalništvo in informatiko, Tržaška 25, Ljubljana, Slovenija E-pošta: marko.bratina@sava.si, andrej.dobnikar@fri.uni-lj.si, uros.lotric@fri.uni-lj.si Povzetek. V prispevku so analizirane možnosti uporabe mer, ki izhajajo iz teorije informacij pri iskanju znacilk in modeliranju dinamičnih sistemov iz Časovnih vrst. Pri iskanju znacilk smo uporabili informacijsko-teoreticno zasnovano analizo neodvisnih osi in metodo najznacilnejših projekcij. Pri gradnji modela smo v nevronski mrezi, dvoplastnem perceptronu, kot kriterijsko funkcijo poleg najpogosteje uporabljane minimizacije povprecne kvadratne napake uporabili še maksimizacijo informacijskega potenciala. Rezultati napovedovanja na vec casovnih vrstah so pokazali, da se predlagane metode iskanja znacilk v vecini primerov obnašajo bolje kot klasicšne metode, medtem ko informacijski potencial kot kriterijska funkcija upocšasni konvergenco. Ključne besede: modeliranje casovnih vrst, predprocesiranje, nevronske mreze, teorija informacij Time-series modeling using information-theory techniques Extended abstract. The paper analyzes the possibility of using measures originating from the information theory in modeling dynamic systems from time series. Our modeling was based on a multilayered perceptron neural network into which several information measures were integrated. The presented information measures are based on the Shannon's definitions of entropy and the divergence given in Eqs. (1) and (2). To simplify the extremely computationally intensive methods, the Renyi's counterparts defined by Eqs. (6) and (7) are preferred instead. For example, the information potential given in (9) can be calculated even without integrations. Besides the information potential, the Renyi's approximations of the measures presented by Eqs. (3) and (5) are very useful in modeling. The measures originating from the information theory were included in two ways: in the data preprocessing and as the criterion function during the process of neural network learning. The basic aim of preprocessing is to reduce the number of inputs to the model and/or reformulate them and thus improve the model performance and generalization ability. Besides classical methods for feature extraction, including the appropriate number of the last values in a time series that results in the best model (H1), the greedy algorithm for finding the combination of inputs, enabling the best model performance (H2) or making a linear combination of original inputs using the principal component analysis (PCA), two methods based on the information theory were tested: the independent component analysis (ICA) and the maximally discriminative projection (MDP). The neural network learning is usually based on minimizing the mean square error E (E). As an alternative, maximizing information potential VR(E) is considered as the learning criterion. The proposed ideas were tested on five different time series. Fig. 1 gives results for predicting the future value in the time series, while average results are summarized in Table 1. The model performance is estimated by several measures: the nor- malized root mean square error (NRMSE), the normalized information potential (NIP) and the number of model parameters (Npar). The ICA preprocessing combined with learning based on information potential yields results comparable to method H2. The latter is, however, computationally extremely intensive. In Fig. 2 and Table 2, classification of the predicted values of time series is presented. As a performance measure, the proportion of correct classification (Pok) is added. The ICA and MDP methods perform equally well in this case. We show that the methods based on the information theory can be efficiently used in retrieval of relevant features from data. Besides, the usage of the information potential as a criterion function in the neural-network learning process is also promising. Key words: time-series modeling, pre-processing, neural networks, information theory 1 Uvod Splošni modeli dinamičnih sistemov poskušajo izluščiti pomembne lastnosti procesov neposredno iz obnašanja merljivih parametrov skozi Cas. Temeljijo na predpostavki, da bodo značilnosti, opaZene v preteklosti, obstajale še naprej. Pred samim modeliranjem vhodne podatke pogosto pripravimo, na primer s filtriranjem ali metodami za iskanje znacilk, s Cimer poskušamo poenostaviti modeliranje. Metode predprocesiranja so zelo različne - od popolnoma statističnih do takšnih, ki temeljijo na uspešnosti končnega modela [6]. Večina splosno uveljavljenih metod za analizo in modeliranje sistemov v veliki meri predvideva statistiko drugega reda. V zadnjem času se za reševanje čedalje zahtevnejsih problemov uveljavljajo rešitve, ki presegajo njene omejitve [1]. Zmogljivosti sodobnih računalniških sistemov so omogočile uporabo idej s področja informacijske teorije tudi pri modeliranju dinamičnih sistemov [2]. Osnove teorije informačij izhajajo s konča prve poloviče 20. stoletja, ko je Shannon postavil matematičšno teorijo za obravnavanje temeljnih vidikov komunikačijskih sistemov. Definičije osnovnih mer, to je entropije, divergenče in povprečšne medsebojne informačije, izhajajo iz verjetnostne teorije in statistike. Vsaka mera v svojem kontekstu opredeljujejo količšino informačije v naključnih spremenljivkah, zato jih lahko učinkovito uporabimo tudi na drugih področjih. Na primer, entropijo in divergenčo kot čenilno funkčijo pri modeliranju s splošnimi modeli [3], povprečno medsebojno informačijo pa kot mero za medsebojno povezanost podatkov [4, 5]. V delu se bomo osredinili na problem izbiranja vplivnih podatkov za potrebe modeliranja in na problem samega modeliranja nelinearnih dinamičšnih sistemov z metodami, zasnovanimi na teoriji informačij. V drugem poglavju bomo predstavili osnovne končepte in razsširitve informačijske teorije, v tretjem poglavju pa bomo nakazali, kako jih je mogoče uporabiti pri pred-pročesiranju podatkov. V četrtem poglavju bomo na kratko predstavili, kako lahko v dvoplastnem perčep-tronu kot kriterijsko funkčijo uporabimo mero, ki izhaja iz teorije informačije. Predstavljeni končepti bodo nato v petem poglavju tudi praktičšno ovrednoteni na problemih napovedovanja znanih čšasovnih vrst. Nazadnje bodo v za-ključšku povzete prednosti in slabosti uporabe končeptov informačijske teorije v primerjavi z ze uveljavljenimi. 2 Mere informacijske teorije Osnovni meri, ki izhajata iz teorije informačij, sta entropija in divergenča [2]. Shannonova entropija, ki meri nedoločenost naključnega vektorja X z gostoto verjetnostne porazdelitvep(x), je definirana kot N H(X) = - Jp(x) logp(x)dx . (1) Podobno Kullbačk-Lieblerjeva divergenča meri podobnost med pravo gostoto porazdelitve p(x) naključne spremenljivke X in njeno očeno r(x), ,p(x) r(x) /p( x) p(x) log z tdx . r(x) (2) N J (X) = £ p(x);n pi(x') Y,H(Xi) - H(X) (3) i=1 Kadar nas zanima količšina informačije, ki jo o naključšni spremenljivki Y z verjetnostno porazdelitvijo q(y) vsebuje naključni vektor X, je ta mera kar enaka medsebojni informačiji I(X; Y) = H(X) + H(Y) - H(X,Y) (4) = H(X) - H(X|Y) , (5) pri čemer je H(X, Y) skupna nedoločenost vektorja X in spremenljivke Y, H(X|Y) pa povprečšna nedoločšenost vektorja X, če poznamo spremenljivko Y. Računanje zgornjih mer je zaradi integralov, ki nastopajo v definičijah, izredno zahtevno. Zato se namesto Shannonove entropije pogosto uporablja Reny-ijeva entropija [7], Hr(X) = 1 — a log / pa (x)dx (6) Medsebojno povezanost spremenljivk Xi z gostotami verjetnostnih porazdelitev pi(xi), ki sestavljajo vektor X = (X1,..., XN)T , lahko očenimo z mero [7] ki ji je v limiti a ^ 1 enaka. Analogno v limiti a ^ 1 Renyijeva divergenča DR(p; r) = T kg |p(x) (p|y) 0 1 dx (7) preide v Shannonovo. V nadaljevanju se bomo omejili na obravnavanje kvadratnih Renyijevih mer s parametrom a = 2. V nasprotju s Shannonovimi merami nobena od relačij (4) in (5) ne velja za Renyijeve mere [5], zato so v nadaljevanju mere, ki izhajajo iz omenjenih relačij, označšene z apostrofom. V praksi se gostota verjetnostne porazdelitve zvezne spremenljivke največškrat očeni z metodo Parzenovega okna [8] K p(X) = K-1Y] G (x - xfc) , (8) k= 1 kjer so xi, i = 1,..., K naključne vrednosti (meritve) vektorja X in G(x) izbrana jedrna funkčija. Pogosto za jedrno funkčijo uporabimo Gaussovo funkčijo G(x) = 1^=1 Ghi (xi), kjer je Ghi (xi) enodimenzionalna Gaussova porazdelitev. Ena od moznosti je, da širino porazdelitve določimo s Silvermanovo očeno [9] hi = 1,06aiK-0'2, kjer je ai standardno odstopanje spremenljivke Xi od povprečja. Ko v definičijo Renyijeve entropije (6) vstavimo nastavek Parzenovega okna (8) ter uposštevamo, da je integral produkta dveh Gaussovih funkčij Gaussova funkčija razlike srednjih vrednosti osnovnih Gaussovih funkčij z dvojno variančo, vidimo, da Renyijev informačijski po-tenčial Vr(X) = p2 (x)dx K Gh(x - xi^ dx i=1 K 1 2 GhV2 (xj - Xi) (9) in s tem entropijo HR(X) = — log VR(X) izračunamo brez računsko poZrešnih integracij [7]. 3 Iskanje značilk Z metodami za iskanje znacilk poskušamo iz vhodnih podatkov izluščiti najpomembnejse značilnosti ter s tem poenostaviti modeliranje in izboljšati odzivnost modelov in njihovo sposobnost posplosevanja. Metode se med seboj ločijo po kriterijski funkciji, s katero izbiramo značilke - ta je lahko zasnovana na lastnostih značilk ali pa kar na uspesšnosti modela. 3.1 Klasične metode Značilne spremenljivke lahko izbiramo s pomočjo spektralne in kovariančne analize, ali pa z bolj ali manj intenzivnim preiskovanjem mogočih naborov z različnimi hevrističšnimi metodami ali iskalnimi postopki, na primer evolučijskimi [10]. V nadaljevanju se bomo omejili na dve hevristični metodi. Pri prvi (H1) bomo za značilke izbrali določeno število zadnjih spremenljivk v časovni vrsti. Število bomo določili tako, da bo modeliranje čim uspesnejse. Pri drugi metodi (H2) bomo nabor značilk gradili postopoma. V vsakem koraku bomo med značšilke vključili tisto od preostalih spremenljivk, ki bo skupaj z ze izbranimi značilkami pripeljala do najboljšega modela. Značšilke lahko sestavimo tudi kot linearne kom-binačije vhodnih spremenljivk. To omogoča metoda glavnih osi (ang. Prinčipal Component Analysis, PCA). Gre za matematični postopek [11], v katerem se iz osnovnih vhodnih spremenljivk sestavi manjše število neodvisnih značilk ali glavnih osi. Prva glavna os je postavljena tako, da ima največjo mogočo variančo. Podobno je vsaka naslednja glavna os postavljena tako, da ima največjo variančo na preostalih podatkih. Ponavadi pri nadaljnjem modeliranju uporabimo nekaj prvih glavnih osi. V nadaljevanju smo izbrali glavne osi, pri katerih je varianča večja od 1 % varianče glavne osi. 3.2 Analiza neodvisnih komponent Analiza neodvisnih komponent (ang. Independent Component Analysis, ICA) [12] predvideva, da lahko izmerjene signale zk zapišemo kot linearno mešaničo statistično neodvisnih signalov, sk, zk = Ask. Z analizo neodvisnih komponent zšelimo poiskati tako matriko B, da bodo značilke xk = Bzk = BAsk kar najboljsi priblizek statistično neodvisnih signalov sk. V prvem koraku s transformačijo z'fc = D-1/2VTzk, v kateri je D diagonalna matrika lastnih vrednosti, V pa matrika lastnih vektorjev avtokorelačijske matrike vektorjev zk, poskrbimo, da je avtokorelačijska matrika vektorjev z'fc identiteta. V drugem koraku nato poiščemo rotacijsko matriko R, ki nam da tak vektor signalov xk = Rzj., pri katerem je izbrana kriterijska funkcija optimalna. (Če želimo, da imajo vektorji xk manjše število komponent kot vektorji zk, na slednjih pred opisanim postopkom uporabimo metodo glavnih osi. Obstaja več metod analize neodvisnih komponent [13, 14], nekaj jih temelji na teoriji informacije. (Če vektorje xk razumemo kot točke naključnega vektorja X, potem lahko za kriterijsko funkcijo uporabimo Renyijevo medsebojno informacijo vektorja, ki izhaja iz (7). Ker pa je njen izracun casovno zahteven, se ponavadi uporabi kar zveza (3), v kateri namesto Shannonove uporabimo Renyijevo entropijo [5]. Skupna entropija vektorja X je invariantna na rotacije [15] , zato lahko zadnji clen v (3) izpustimo in minimiziramo samo vsoto entropij posamicnih signalov JR(X) = ^HR(X®). Rotacijsko matriko, ki minimizira mero J'R(X), išcemo iterativno, na primer po postopku najhitrejšega sestopa, po katerem je sprememba rotacijske matrike enaka AR = -n J (X)/d R [5]. 3.3 Metoda najznačilnejših projekcij Tako kot pri analizi neodvisnih komponent tudi pri metodi najznačilnejših projekčij (ang. Maximally Disčrimina-tive Proječtions, MDP) isčemo značilke xk, ki so linearne kombinačije osnovnih meritev zk z dvofaznim postopkom [16]. Bistvena razlika med metodama je v zasnovi kriterijske funkčije. Medtem ko je pri analizi neodvisnih komponent ta odvisna samo od vhodnih podatkov, jo pri metodi najznačilnejših projekčij določajo odvisnosti med vhodnimi podatki in ustreznimi izhodi iz modela. Opisanemu končeptu ustreza maksimizačija mere medsebojne informačije med vhodnimi in izhodnimi spremenljivkami. Računanje prave Renyijeve medsebojne informačije je preveč zahtevno, zato se uporablja priblizek, ki izhaja iz (5), v kateri Šhannonove entropije zamenjamo z Reny-ijevimi. Ta mera se se dodatno poenostavi v primeru, ko zelimo z modelom vhodne podatke uvrščati v C vnaprej podanih razredov. Takrat lahko zapišemo Ir (X; Y) = Hr (X) — Hr(X|Y ) = Hr(X) —V ^Hr(X|Y= c) , (10) i—' n c=1 kjer je nc število podatkov, ki se uvršcajo v razred c, HR(X| Y = c) pa nedolocenost vektorja X pri uvršcanju v ta razred. Pri racunanju entropije HR(X) uporabimo vse vzorce, pri racunanju entropije HR(X|Y = c) pa le vzorce, ki spadajo v razred c. Tako kot pri analizi neodvisnih komponent rotacijsko matriko popravljamo iterativno glede na vrednosti gradienta kriterijske funkcije [16]. 4 Modeliranje z nevronskimi mrežami Nevronske mreže spadajo med splošne nelinearne modele, ki izbrano vrednost v časovni vrsti povezujejo s predhodnimi vrednostmi. Nevronska mreža dvoplastni perceptron [17], ki smo jo uporabili pri analizi, ima v skriti plasti MH nelinearnih nevronov, v izhodni plasti pa MO linearnih nevronov. Odziv modela lahko opišemo z enačbama xf = tanh (Wh xfc + bH) in (11) xO = Wo xH + bo . (12) Z učenjem na parih vhodno-izhodnih vzorcev (xk, dk), k = 1,..., K proste parametre modela, utezi WH in WO ter pragove bH in bO nastavimo tako, da optimiziramo izbrano kriterijsko funkcijo na podlagi napak ek = dk — xO, k = 1,..., K med dejanskimi Ponavadi minimiziramo 1_ V^K T No K ek ek. in izracšunanimi vrednostmi. povprečno kvadratno napako £(E) Med merami, ki izhajajo iz informacijske teorije, je primerna minimizacija nedolocšenosti napake. Ob predpostavki, da napake ek tvorijo nakljucni vektor E, se za kriterijsko funkcijo najvecškrat uporablja minimizacija entropije oziroma maksimizacija informacijskega potenciala VR(E), podanega v (9). Pri minimizaciji povprecšne kvadratne napake z gra-dientnimi metodami gradiente izracšunamo po znanem vzvratnem postopku [17]. Po analognem postopku je mogocše izracšunati tudi gradiente informacijskega potenciala [3]. Informacijski potencial ni odvisen od povprecja porazdelitve napak, zato se lahko zgodi, da po koncani optimizaciji povprecje napak ne bo nic. Ker pa so izhodni nevroni linearni, lahko anomalijo odpravimo tako, da po koncanem ucenju ustrezno nastavimo pragove izhodnih nevronov [7]. 5 Eksperimentalno delo Omejili se bomo na modeliranje diskretnih cšasovnih vrst, vzorcšenih v enakomernih cšasovnih presledkih. Za napovedovanje smo izbrali pet casovnih vrst: povprecno letno število soncevih peg (Pege), kaoticno logisticno preslikavo (LP), kaoticno casovno vrsto Mackey-Glass (MG) ter tecaj delnice podjetja Sava (Sava) in vrednosti borznega indeksa SBI20 (SBI) v obdobju od zacetka aprila 2007 do konca junija 2008. Modelirali smo na dva nacšina: pri prvem smo napovedovali vrednost cšasovne vrste v naslednjem koraku, pri drugem pa smo pricakovane spremembe vrednosti v casovni vrsti uvršcali v pet razredov. V vsaki casovni vrsti je bilo 308 podatkov, od katerih smo jih prvih 80 % uporabili za modeliranje, zadnjih 20 % pa za testiranje in primerjavo uspesnosti modelov. Da smo omejili velikost modelov, so le-ti smeli za napoved naslednje vrednosti uporabiti 12 predhodnih vrednosti, poleg tega pa sštevilo prostih parametrov modela ni smelo presegati 40 % števila podatkov v casovnih vrstah. Da bi se izognili naključnim zaustavitvam optimizacije v lokalnih minimumih, smo vsak model zgradili desetkrat, vsakič z naključno določenimi začetnimi parametri. V nadaljevanju so predstavljeni najboljši rezultati na testnih mnozičah. 5.1 Napovedovanje vrednosti Za napovedovanje naslednje vrednosti v časovni vrsti z dvoplastnim perčeptronom so bili uporabljeni štirje različni načini izbiranja značilk in dve kriterijski funkčiji. Uspešnost napovedovanja je pri vsakem modelu podana s tremi čenilnimi funkčijami: s koren-jeno povprečno kvadratno napako, normalizirano na standardno odstopanje časovne vrste od povprečja (ang. Normalized Root Mean Squared Error), NRMSE = \JE(E)/ct, z normaliziranim informačijskim poten-čialom (ang. Normalized Information Potential), NIP = VR(E)/max{VR(E)} in s stevilom prostih parametrov modela (Npar). Rezultati so v obliki grafikonov predstavljeni na sliki 1. Na grafikonih vidimo, da sta vrednosti w0.6 C/5 §0.4 ^ 0.2 0 10.9 0.8 60 J 40 t? 20 0 il L tu Ulj M lllli nil Lhn 111 lil J" Pege LP MG Sava SBI ■ HI hH2 IPCA i=i ICA