- поточний (оперативний) прогноз, що орієнтований на строк упередження, в межах якого можливі кількісні зміни параметрів об'єкта дослідження;

- довготерміновий (довгостроковий) прогноз - в межах його можливі якісні (структурні) зміни об'єкта дослідження та оточуючого середовища, в якому він функціонує;

- середньотерміновий (середньостроковий) прогноз - проміжний між коротко- та довготерміновим з можливою перевагою очікуваних кількісних змін над якісними;

- наддовготерміновий (наддовгостроковий) прогноз, в межах якого можливі вельми суттєві якісні (структурні) зміни досліджуваного об'єкта та середовища його функціонування, через що доцільно говорити лише про найзагальніші перспективи прогнозу.

Безпосередньо для розробника прогнозу вельми актуальне питання ступеня формалізації методу прогнозування. За цією ознакою всі методи прогнозування поділяються на дві категорії: евристичні та математичні, в основі кожної з яких лежить один з двох альтернативних підходів до розв'язку задач прогнозування.

Евристичні методи (від грецького слова "heurёka" - робити відкриття, або від більш відомої, завдяки історичному анекдотові, репліки Архімеда “Еврика!” - знайшов) започатковані ще в Древній Греції та Римі. Ці методи не потребують формального визначення проблеми і, відповідно, побудови алгоритму, що дає строгий й цілком інтерпретований розв'язок. Визначальна риса евристичного прогнозування - активне використання в тій чи іншій формах знань, досвіду та інтуїції експертів - висококваліфікованих фахівців з відповідної предметної галузі. Отриманий евристичний прогноз завжди інтегрований з конкретною особою - експертом чи групою експертів і має суб'єктивний характер.

В свою чергу евристичні методи поділяються на експертні (інтуїтивні) та аналітичні.

Експертні методи ґрунтуються виключно на оцінках експертів, зроблених відносно проблеми, яку вивчають. При цьому механізм продукування цих оцінок лишається невизначеним. Як правило, він невідомий навіть самому експерту, носить виключно індивідуальний, особистий характер і не може бути повторений чи відтворений кимсь іншим.

Натомість аналітичні методи прогнозування базуються на логічному (теоретичному та емпіричному) аналізі усього наявного масиву відомостей про об'єкт прогнозування, наслідком чого є побудова певної логічної схеми (моделі). Однак треба мати на увазі, що формування цієї логічної схеми звичайно проходить в умовах неповної визначеності щодо властивостей об'єкту прогнозування через завжди існуючий брак вхідної інформації. Недостатня інформація домислюється експертом відповідно до своїх вподобань та пріоритетів, тому отримана схема є суб'єктивною і, хоч дозволяє побудувати відповідну технологію прогнозування, придатну для багаторазового використання, отримані нею результати треба завжди сприймати як індивідуальний прогноз, зроблений виключно через призму особистості конкретного експерта.

Тому головна проблема евристичних методів - це визначення і врахування під час обробки експертних оцінок рівня компетентності експерта, його здатності робити достовірні припущення (прогнози) щодо конкретних показників та перспектив розвитку об'єкту (процесу) прогнозування.

Математичне прогнозування базується на використанні математичної моделі об'єкту (процесу) прогнозування, побудованої за ретроспективними даними (тобто, отриманими до певного моменту часу в минулому), для обчислення його характеристик (станів) в довільний заданий момент.

Аргументами математичної моделі можуть бути контрольовані фактори (інформаційні ознаки), зміни кількісних значень яких впливають на значення вихідної змінної. Якщо до складу факторів входить час, об'єкт прогнозування належить до класу динамічних об'єктів, якщо тільки час - до часових рядів.

Сам процес математичного прогнозування можна умовно поділити на чотири головних етапи:

збирання та підготовка вихідних даних;

пошук та обґрунтування структури математичної моделі об'єкту (процесу) прогнозування (задача структурної ідентифікації моделі);

оцінювання параметрів моделі за вихідними даними (задача параметричної ідентифікації моделі);

перевірка якості моделі та прогнозування.

Ключовим етапом у наведеній процедурі є другий, зміст якого треба було б доповнити одним обов'язковим припущенням, невиконання якого лишає прогноз будь-якого сенсу: модель процесу (об'єкту) побудована за результатами аналізу ретроспективних даних, не змінюється впродовж часу упередження, тобто з поточного моменту до моменту, для якого обчислюють прогноз. Якщо це припущення виконується, точність прогнозу цілком залежить від рівня адекватності моделі об'єкту прогнозування, тобто від якості розв'язання задачі структурно-параметричної ідентифікації моделі.

Через те, що для моделі неминучі певні спрощення порівняно з оригіналом (об'єктом моделювання), точність прогнозу є принципово обмеженою, у зв'язку з чим постає питання про можливий (чи доцільний) інтервал упередження (прогнозу). Як вже зазначалося, у більшості прикладних застосувань розрізняють довго-, середньо- чи короткострокові прогнози. Для останнього звичайно (якщо це можливо в принципі) використовують математичне моделювання. Евристичні методи використовують для довгострокових прогнозів. Крім того, їх можна використовувати у разі якісної форми опису об'єкта прогнозування, тоді як математичні для побудови кращої моделі потребують кількісної форми подання вихідних даних.

В цьому випадку як вхідні дані, так і обраховані прогнозні (модельні, вихідні) значення мають кількісну форму представлення. Тому група відповідних методів, що забезпечують отримання кількісного прогнозу, дістала назву методів кількісного прогнозування. Типовими представниками цієї групи є методи прогнозування часових рядів - впорядкованої у часі множини дискретних відліків певного параметру, величина якого змінюється у часі. Метою прогнозування часових рядів є визначення майбутньої поведінки контрольованого параметру на певному інтервалі упередження, розташованому справа від відмітки поточного моменту на часовій координаті. У наступних розділах буде розглянуто ряд методів прогнозування часового ряду, які можна вважати традиційними для цього класу задач.

Як зазначалося вище, традиційний підхід до кількісного математичного прогнозування базується на умові однаковості, незмінності математичної моделі (так звана "жорстка" модель [20-23]), як на інтервалі аналізу, тобто на ретроданих, так і на інтервалі упередження. Це означає, що фактично розглядається задача екстраполяції ретроспективної тенденції розвитку процесу, представленого сукупністю елементів часового ряду, на майбутнє, що принципово виключає можливість прогнозу змін, "стрибків", будь-яких інших кардинальних відступів від сформованої за ретроданими "жорсткої" моделі ряду. Тому для реалізації середньо- та довгострокових прогнозів, а також прогнозів, що потребують якісної форми опису та представлення прогнозної інформації, потрібне застосування нових методологічних принципів та підходів. Спроби їх пошуку та застосування достатньо детально описані в літературі [16,17].

6.2 Часові ряди. Структура, моделі

Однією з ознак, за якою здійснюється класифікація моделей, є врахування у моделі фактору часу. За цією ознакою моделі поділяються на статичні та динамічні. Перші описують лише ті закономірності та співвідношення між параметрами ОД, існування яких не залежить від часу (не змінюється у часі). Динамічні моделі дозволяють показати процес функціонування та розвитку ОД, аналізувати тенденції й закономірності цього розвитку, дослідити вплив сукупності різних факторів на механізм цього розвитку у часі. Проте процедура ідентифікації динамічних моделей висуває специфічні вимоги щодо формату подання вихідних даних: вони мають бути представлені у формі часових рядів.

Часовий або динамічний ряд - упорядкована у часі послідовність дійсних чисел , що характеризує поведінку ОД впродовж терміну його спостереження. Чисельні значення кожного окремого елементу часового ряду {}, називають його рівнями. Їх вимірюють у послідовні моменти часу (звичайно через рівні проміжки, так званий регулярний ряд), через що порядок розташування елементів часового ряду дуже суттєвий.

Часовий ряд - особливий тип даних, який містить інформацію про плинність значень станів або характеристик ОД у часі. Ряд послідовних елементів може відображати деяку закономірність, тенденцію у функціонуванні чи розвитку ОД, наприклад, попередні значення ряду впливатимуть на зміну поточних чи майбутніх його значень. Таким чином, припускається наявність зв'язків або залежності між елементами одного й того ж ряду, яку, зокрема, можна зафіксувати у формі ММ часового ряду. Якщо вірогідне збереження цих закономірностей (зв'язків, залежностей) в майбутньому, реальною є можливість на основі побудованої за ретроданими (тобто даними, що описують поведінку часового ряду в минулому) моделі часового ряду передбачити його значення на певний період упередження.

Звичайно при дослідженні часових рядів в їх структурі намагаються окремо виділити детерміновану та стохастичну складові, процедури моделювання та прогнозування яких принципово відмінні. Найчастіше при цьому спираються на адитивну модель представлення часового ряду:

, . (39)

де стохастична складова часто традиційно асоціюється із випадковою похибкою Е , у відліках якої відсутня автокореляція:

а математичне очікування якої .

В свою чергу в детермінованій складовій традиційно виділяють три адитивних компонента: тренд , сезонну та циклічну , тобто

. (40)

Ці компоненти не мають єдиних загальноприйнятих визначень й змістовно залежать від конкретної предметної сфери застосування (іноді трендом називають загалом всю детерміновану складову моделі). В соціальних дисциплінах названі компоненти найчастіше ідентифікуються наступним чином.

Тренд - базовий (основний) детермінований компонент часового ряду, який характеризує вплив довгострокових постійно діючих факторів, повільний та інерційний.

Сезонний компонент відображає достатньо поширену як серед природних, так і штучних процесів регулярну повторюваність явищ: метеорологічних, кліматичних, соціально-економічних, демографічних й т.п., які не мають точної періодичності та сталості амплітуд. Приклади: зміни кожноденної завантаженості міського пасажирського транспорту, коливання трафіку (кількість з'єднань) у телефонній мережі протягом доби, сезонна регулярність термінів проведення сільськогосподарських робіт ін.

Циклічний компонент описує нерегулярні, відносно довгі за часом періоди зростання та спаду рівнів часового ряду: зміни рівня радіації Чорнобиля, зміни чисельності популяції під час та після епідемії, зміни показників економіки внаслідок соціально-економічних спадів. Циклічний компонент важко ідентифікувати лише формальними засобами, необхідно залучати додаткову інформацію.

У аналізі часових рядів головну увагу звичайно приділяють виділенню трендового компонента, вважаючи його основою детермінованої складової часового ряду, на яку накладаються інші структурні компоненти моделей (39), (40). Тому далі в цьому розділі будуть розглядатися питання, пов'язані виключно з ідентифікацією трендової моделі часового ряду.

Залежно від сфери прикладного застосування накопичено певний досвід використання математичних моделей та визначено певну сукупність найбільш поширених базових трендових моделей. Зокрема, це показникові, степеневі, обернені та інші моделі [16]:

, , , , ,

, ,... .

Однак, якщо типової базової форми моделі заздалегідь визначити не вдається, найбільш гнучкою з точки зору її апроксимативних можливостей і наступного використання є поліноміальна модель виду:

, (41)

де фактором, що впливає на залежну змінну , є час , який входить до моделі в різних степенях у вигляді відповідних регресорів: . В поліноміальній моделі звичайно використовується дискретний час t, значення якого співпадають з натуральним рядом чисел, тобто цей час слід вважати фіктивною змінною, що в моделях задається безвідносно до реального обліку часової змінної.

7. Прогноз трендового компонента часового ряду

7.1 Поліноміальні трендові моделі. Оцінювання якості, критерій Дарбіна-Уотсона

Поліноміальна трендова модель (41) фактично є рівнянням множинної регресії (11), тому для її ідентифікації цілком застосовні розглянуті в першій частині цього видання методи та процедури регресійного аналізу. Зокрема розширена матриця даних для моделі (41) матиме вигляд:

, (42)

а вектор МНК-оцінок параметрів визначатиметься співвідношенням:

. (43)

Проте на етапах вибору варіанту ММ та її верифікації виникає, завдяки специфічним особливостям часових рядів, можливість використання додаткових способів оцінювання якості варіантів ММ. Найчастіше це стосується аналізу нев'язок (залишків, похибок) , трендової моделі, де - оцінка значення t-ого рівня часового ряду, розрахована за поліноміальною трендовою моделлю (41).

Зокрема, зважаючи на характер обмежень, що накладаються на властивості випадкової складової et рівня часового ряду, а саме, на припущення щодо некорельованості її значень:

видається доцільним перевірити, як виконується це обмеження для нев'язок , що за певних обставин сприймаються у якості оцінок відповідних значень випадкової складової Е. Найбільш відомим і поширеним тестом перевірки моделі на наявність автокореляції між нев'язками є критерій Дарбіна-Уотсона. Критерій застосовується у тих випадках, коли тренд визначається для всього ряду [1,26,27].

Змістовну суть критерію можна пояснити наступним чином. У разі правильної специфікації трендової моделі за даними всього часового ряду, тобто коли структура функції тренду адекватно відображає динаміку змін значень рівнів часового ряду на інтервалі його спостереження, послідовність нев'язок дійсно слід вважати оцінками значень випадкової складової Е. Очевидно, що в цьому випадку автокорельованість нев'язок має бути мінімальною. У протилежному випадку, коли отримана трендова модель є недостатньо точною, в нев'язках присутня, крім випадкової складової , частка детермінованої складової часового ряду [1], яка має високий рівень автокореляції. Тому, наприклад, серія позитивних або негативних знаків у послідовності нев'язок є своєрідним індикатором того, що побудована модель неадекватна вихідним даним. Більш вагомим статистичним підтвердженням цього факту є існування автокореляцій у серії нев'язок, одна з основних причин появи якої - наявність в них залишкового тренду.

Значення статистики Дарбіна-Уотсона розраховується за формулою:

. (44)

З аналізу формули (44) витікає, що якщо між значеннями присутня сильна позитивна автокореляція , де - вибірковий коефіцієнт автокореляції між двома послідовними нев'язками часового ряду, то величина d?0, якщо маємо сильну негативну автокореляцію , то d?4. При відсутності автокореляції , статистика d?2. Таким чином, значення d-статистики Дарбіна-Уотсона знаходяться у межах від 0 до 4 (). На рис. 3 зображено схему розташування граничних значень d-статистики:

Рис.3. Зони розташування граничних значень d-статистики

Крім формули (44) для приблизного розрахунку статистики Дарбіна-Уотсона можна використати простішу формулу:

. (45)

Співвідношення (45) дає найбільш прозору інтерпретацію змісту критерію Дарбіна-Уотсона та його кількісним характеристикам.

Практичне використання критерію Дарбіна-Уотсона засновано на порівнянні величини d, розрахованої за формулою (44), з теоретичними значеннями d1 та d2, які показані на рис.3:

1)при - приймається гіпотеза про відсутність автокореляції між значеннями випадкової величини;

2)при - виявляється позитивна автокореляція;

3)при - виявляється негативна автокореляція;

4)при або - немає достатньої інформації для прийняття рішень, тобто отримуємо зону невизначеності.

Теоретичні граничні значення d1 та d2 для рівня значущості 5% і різних степенів поліному p представлено у табл. 2 (більш повні дані наведено у додатках в [3,4]).

Таблиця 2. Граничні значення d1 та d2 ( р - степінь поліному, n - об'єм вибірки)

n	p=1	p=2	p=3	p=4
	d1	d2	d1	d2	d1	d2	d1	d2
30	1,35	1,49	1,28	1,57	1,21	1,65	1,14	1,74
40	1,44	1,54	1,39	1,6	1,34	1,66	1,29	1,72

Для визначення за відомою (вже ідентифікованою) поліноміальною моделлю виду (41) прогнозних значень тренду достатньо ввести до цієї моделі замість змінної t її кількісне значення l, в якому враховано інтервал упередження L. Так, якщо цей інтервал дорівнює одному, двом, трьом,... відлікам (крокам), тобто L=1,2,3,..., то значення змінної t кількісно дорівнюватиме l= n + L.

В матричній формі оцінка прогнозу тренду визначається виразом

, (46)

де - вектор значень регресорів моделі (41) у точці прогнозу.

Помилка отриманої оцінки трендового прогнозу матиме (за умови адекватності структури моделі вихідним даним та обчисленню коефіцієнтів поліноміальної моделі методом найменших квадратів) дисперсію

. (47)

Слід зазначити, що процедура виділення тренду суттєво ускладнюється в разі наявності в даних часового ряду циклічного компоненту (особливо в коротких чи середніх за обсягом рядах). Зокрема, помилкове включення рівнів, у складі яких присутній циклічний компонент, до сукупності даних, що використовується для ідентифікації трендової моделі, призводить до створення трендової залежності і робить неможливим її використання для прогнозу тренду часового ряду. В цьому випадку рекомендується [2] проаналізувати одночасну сукупну зміну у часі трендового та циклічного компонентів, відділивши від них шумову складову , шляхом згладжування даних часового ряду.

7.2 Згладжування часових рядів зваженим ковзким середнім. Критерій Аббе

Мета згладжування часового ряду - виділення його детермінованої складової , . Одним з найбільш поширених способів згладжування часового ряду є його обробка із застосуванням процедури зваженого ковзкого середнього. В загальному вигляді процедуру обчислення згладженої оцінки для t-ого рівня ряду описує вираз:

, (48)

де , - вагові коефіцієнти згладжуючого фільтру - рівні ряду, що утворюють так званий ковзкий інтервал (КІ) довжиною 2m+1, який отримав назву завдяки своєму покроковому переміщенню вздовж часового ряду у ході виконання згладжування останнього. Покрокове переміщення КІ здійснюється шляхом відкидання його крайнього лівого елементу та приписування до КІ справа найближчого сусіднього рівня часового ряду. При такому зміщенні КІ черговий згладжуваний рівень ряду опиняється у середині КІ, тобто оцінка завжди обчислюється для середини КІ. Через це m перших та m останніх рівнів часового ряду залишається незгладженими. Це явище отримала назву крайового ефекту.

Для отримання згладженої оцінки за даними КІ будується поліноміальна модель виду (41), яка має апроксимувати фрагмент часового ряду в межах ковзкого інтервалу. Для цієї моделі розраховується вектор МНК-оцінок параметрів , який далі використовується для обрахування згладженої оцінки в середній точці КІ. Отримана оцінка є шуканою згладженою оцінкою t-ого рівня часового ряду. За аналогією з (46) маємо:

, (49)

де - середній рядок матриці плану .

Зважаючи, що довжина 2m+1 ковзкого інтервалу та порядок р поліному (41) залишаються незмінними впродовж всієї обробки часового ряду, незмінною буде і матриця плану , тому матричний добуток у правій частині виразу (49), пов'язаний лише з параметрами m та р, можна обрахувати заздалегідь:

. (50)

Зробивши відповідну заміну у співвідношенні (49):

,

отримуємо матричну форму запису виразу (48), тобто - вектор вагових коефіцієнтів згладжуючого фільтру. Властивість сталості значень вагових коефіцієнтів при заданих параметрах m та р фільтрів дозволяє виконати попереднє обрахування ваг для різних сполучень m та р та внести отримані значення векторів до довідкових розділів спеціалізованої літератури [4].

Для практичної реалізації процедури згладжування актуальне питання вибору параметрів m, р згладжуючого фільтру. Звичайно воно розв'язується шляхом підбору значень m, р за допомогою критерію якості згладжування. Методичною базою цього перебору може бути викладена у розділі 7.1 ідея перевірки автокорельованості нев'язок , . Саме ця ідея реалізована у критерії Аббе, який можна застосувати для вибору параметрів фільтру:

, (51)

де - середнє арифметичне значення нев'язок.

Якщо , то приймається гіпотеза некорельованості нев'язок, тобто якісного згладжування. Теоретичне значення для n>60 розраховується за формулою:

, (52)

де - б-квантиль нормованого нормального розподілу. При значення для трьох найбільш застосовних рівнів значущості б наведено в [18].

Зазначимо, що оцінку дисперсії нев'язки, отриману для найближчого до 1 значення статистики Аббе, можна вважати найбільш точною оцінкою дисперсії стохастичної складової еt часового ряду:

. (53)

7.3 Оцінка точності прогнозу за ретроданими

Нажаль розглянуті вище критерії Дарбіна-Уотсона (44) та Аббе (51) дають достатньо чітке уявлення про апроксимативні якості моделей на інтервалі спостереження, але не відповідають на питання про якість прогнозу. Більш-менш об'єктивну відповідь на це питання можна отримати використовуючи ретропрогноз, тобто прогноз, знайдений за трендовою моделлю, побудованою для фрагменту часового ряду , , , де L - інтервал прогнозу (звичайно ), , при цьому було б бажано, щоб для n1 та n2 за можливістю виконувались умови: , n2>10. Розрахованому в цьому випадку прогнозному значенню рівня слід співставити відомий елемент вихідного ряду, визначивши таким чином нев'язку ретропрогнозу, що складатиме , де нижній індекс R означає приналежність до ретропрогнозу. Наступна прогнозна оцінка розраховується за подовженим на одиницю фрагментом ряду , , відповідна нев'язка дорівнюватиме . Подібні підрахунки продовжуватимуться доти, доки на фрагменті , не буде обчислений останній ретропрогноз . Слід наголосити, що при ідентифікації моделей на фрагменті часового ряду структура цих моделей має співпадати із структурою трендової моделі, визначеною за даними всього часового ряду (відповідно до розділу 7.1).

За отриманим рядом нев'язок ретропрогнозу можна проаналізувати їх сталість або наявні тенденції до змін, спробувати оцінити дисперсію прогнозної оцінки . Для цього знайдемо оцінку дисперсії нев'язок ретропрогнозу. В якості її можна прийняти середній квадрат рівнів наведеного вище ряду нев'язок або оцінку, обраховану для цього ряду методом експоненціального згладжування. В останньому випадку, застосовуючи рекурентну процедуру

, , (54)

де в якості початкової оцінки приймається значення , розраховане за формулою (47), отримуємо оцінку дисперсії ретропрогнозу для кінця часового ряду.

Зважаючи на те, що дисперсія ретропрогнозу складається з двох частин, одна з яких - дисперсія стохастичної складової часового ряду , а друга - саме дисперсія прогнозу на інтервалі упередження L, знаходимо:

, (55)

де визначається за формулою (53).

7.4 Метод гармонічних вагових коефіцієнті

Метод гармонічних вагових коефіцієнтів був запропонований польським статистиком З.Хелвигом для прогнозування часових рядів, що не мають сезонних та циклічних коливань. Основна ідея методу - обробка рівнів часового ряду таким чином, що більш пізнім його рівням надається більша вага. Це досягається через введення у процедуру обробки спеціальних вагових коефіцієнтів.

Переваги методу гармонічних ваг у порівнянні з іншими методами, де також використовуються зважування рівнів часового ряду, у тому, що його застосування не потребує ніяких припущень щодо виду тренду. Розглянемо цей метод більш докладно.

Нехай маємо часовий ряд , який можна розкласти на невипадкову функцію від часу (тренд) та стаціонарний випадковий компонент: . Треба визначити оцінку прогнозу рівня для ()-ого елемента часового ряду, де - інтервал упередження. Процедура прогнозу рівня складається з трьох етапів: згладжування вихідного часового ряду , обчислення середньозваженого приросту згладжених значень , на інтервалі, який дорівнює інтервалу упередження , визначення оцінки прогнозу.

Згладжування вихідного часового ряду , зазвичай реалізується методом зваженого ковзкого середнього першого порядку (р=1) з ковзким інтервалом (КІ) довжиною три або п'ять елементів (тобто m=1 або m=2). Для запобігання скороченню згладженого часового ряду на 2m елементів внаслідок так званого крайового ефекту (відкидання m перших та m останніх елементів ряду при застосуванні для обчислення згладжених оцінок формули 48)) використаємо формулу

, (56)

в якій на відміну від формули (48) використовуються вагові коефіцієнти , з подвійною індексацією, що дозволяє, залежно від значення k (), обчислити згладжені оцінки для будь-якого елемента КІ. Використовуючи цю можливість, за формулою (56) для m=1 при , обраховується згладжена оцінка першого рівня часового ряду, при , k =1 - останнього, у випадку m=2 при , - дві перші оцінки та при , - відповідно дві останні. При k=0 обчислюються згладжені оцінки у середині КІ, тобто формула (56) співпадає із (48).

Обчислимо послідовність приростів для пар рівнів, розділених L періодами:

, (57)

У відповідність кожному приросту поставимо гармонічну вагу, яка обчислюється за формулою:

, (58)

або через рекурентну процедуру, що задається співвідношенням:

+, де . (59)

Сума гармонічних ваг дорівнюватиме

. (60)

Вводячи нормуючий множник 1/S(L), розрахуємо сукупність гармонічних коефіцієнтів:

,, (61)

сума яких дорівнюватиме одиниці (так звана умова незміщеності системи коефіцієнтів):

; (62)

Розрахуємо середньозважені прирости для різних значень L:

, (63)

Прогнозні значення рівнів ряду з інтервалом упередження L обчислюються за формулою:

, (64)

Можливі варіанти виразів для обчислення рівнів прогнозів на інтервалах L=1,2,3 наведено нижче у таблиці.

1-й варіант	2-й варіант	3-й варіант

8. Прогнозування стохастичної складової часового ряду. Метод авторегресійних моделей

8.1 Модель авторегресії: основні поняття

При наявності лінійної залежності між випадковими величинами в послідовності зустрічаємося з явищем авторегресії - регресії, пояснюючими змінними якої є попередні рівні динамічного ряду [1,2,3].

Моделлю авторегресії називається модель стаціонарної послідовності, що відображає значення показника у вигляді лінійної комбінації скінченого числа попередніх значень цього показника та адитивної випадкової складової:

. (65)

Лаг (або часовий лаг, лаг запізнення) - проміжок часу , за який зміна аргументу призведе до зміни результативного показника. Наявність запізнення означає, що вплив змінної x на змінну y не проявляється негайно, а розтягується на деякий проміжок часу. Значення р в моделі (65) визначає порядок авторегресії. Якщо в моделі часового ряду використовується дискретний час t, значення якого співпадають з натуральним рядом чисел, максимальний лаг визначає порядок авторегресії р.

Розглянемо, чим може бути корисна модель авторегресії для прогнозу часових рядів. Як зазначалося у п. 6.2, часовий ряд можна розкласти на детерміновану функцію часу , котра частіше за все асоціюється із трендом, та випадкову складову .

Припустимо, що, приймаючи введені в п. 6.2 певні обмеження щодо властивостей складової в елементах часового ряду , , її можна вважати реалізацією дискретного білого шуму, отриманого внаслідок дискретизації "чисто випадкового" стаціонарного стохастичного процесу з нульовим математичним очікуванням. Тоді , де - інформативна складова i-ого елемента часового ряду. Обчислені після виділення тренду нев'язки , в загальному випадку відрізняються від , бо в своєму складі містять ще один додатковий компонент - помилку апроксимації трендовою моделлю інформативної складової часового ряду. Частіше за все можна інтерпретувати як наслідки присутності в часовому ряді певної сукупності мікроциклових компонентів, що не можуть бути описані трендовою моделлю. Тому фрагменти сукупності , являють собою сильнокорельовані послідовності, які за своїми частотно-часовими властивостями значно ближче до рядів , , ніж до нев'язки , з елементами якої їх зв'язує співвідношення

, , (66)

тобто значення нев'язки за певних обставин можна вважати оцінками компоненту .Зважаючи на це, якщо припустити, що послідовність апроксимується моделлю виду

, (67)

і це припущення виявиться слушним, то, підставивши до правої частини моделі відповідні значення (які у даному випадку виступають у якості оцінок рівнів ряду ), можна обрахувати оцінки корегуючої поправки , додавання якої до трендового прогнозу дозволить підвищити точність останнього. Таким чином загальний прогноз часового ряду буде складатися з результатів двох прогнозів:

. (68)

Тобто за допомогою авторегресійної моделі (65) з нев'язок трендової моделі , можна виділити додаткову корисну інформацію щодо апріорно невідомої інформативної складової елементів часового ряду. Однак для цього слід попередньо визначити структуру моделі авторегресії (65), зокрема її порядок р. Відомості про це можна отримати з аналізу автокореляційної функції нев'язки , точніше з вибіркових оцінок значень цієї функції , розрахованих відповідно для дискретних моментів часу k=0,1,2,... .

Для обчислення корегуючої поправки оцінок тренду часового ряду за допомогою авторегресійних моделей необхідно виконати ряд послідовних операцій:

виключити тренд , уточнивши структуру моделі тренду відомими способами;

розрахувати відхилення (нев'язки) моделі тренду;

для часового ряду розрахувати автокореляційні функції ;

знайти порядок авторегресійної моделі та оцінити її параметри ;

використовуючи отримані оцінки параметрів авторегресійної моделі, побудувати корегуючу поправку трендового прогнозу часового ряду;

оцінити ефективність введення корегуючої поправки.

Ці операції можуть повторюватися в процесі уточнення досліджуваних моделей, забезпечуючи за рахунок аналізу значень різноманітних показників властивостей даних та моделей логічні обґрунтування отриманих числових оцінок прогнозу при збереженні припущення про незмінність умов та характеристик часового ряду на період прогнозування.

8.2 Оцінювання структури та коефіцієнтів рівняння авторегресії

Процес авторегресії характеризується тим, що його автокореляційна функція є затухаючою (рис. 4), на відміну від автокореляційної функції білого шуму, яка теоретично має дорівнювати 0 для всіх . При завжди дорівнює 1.

Рис.4. Приклади автокореляційних функцій стаціонарних випадкових процесів



Рис. 5. Фрагмент реалізації низькокорельованого випадкового процесу (а) та оцінка його автокореляційної функції (б)

Графік вибіркових оцінок автокореляційної функції називається корелограмою, яка характеризує рівень лінійного зв'язку між елементами часового ряду, віддаленими один від одного відповідно на 0,1,2,... періоди опитування (дискретизації). Для побудови корелограми необхідно:

1) оцінити дисперсію нев'язки:

; (69)

2) розрахувати вибіркові значення кореляційної функції

, (70)

де n1 задається рівним від третини до половини обсягу сукупності n.

На рис. 5а відображено вихідний ряд білого шуму та на рис. 5б - корелограму шуму. Для гаусовського білого шуму наближено [1,2] можна вказати 95% довірчий інтервал вибіркових оцінок значень : . Його зображено на графіку корелограми (рис. 5б) двома тонкими горизонтальними лініями. Якщо вибіркові оцінки кореляційної функції не виходять за межі вказаного довірчого інтервалу, можна припустити, що дані, за якими розраховані оцінки , є реалізацією дискретного білого шуму.

Одним із важливих етапів побудови авторегресивної моделі є визначення її порядку та обчислення її параметрів. Якщо порядок моделі вибрано не вдало, то процедуру обчислення оцінок можна повторити для моделі авторегресії іншого порядку чи структури. Необґрунтоване підвищення порядку моделі та ускладнення її структури знижує точність оцінок параметрів та якість прогнозу. Водночас недостатня кількість коефіцієнтів моделі й занадто малий порядок авторегресії не дадуть можливість адекватно оцінити динаміку процесу та спрогнозувати його подальші зміни. Тому задача дослідника - вибрати модель авторегресії найменшого порядку за умови забезпечення достатньої точності опису даних та прогнозування.

Для визначення порядку авторегресії використовують автокореляційну функцію . Визначають значущі значення вибіркових оцінок автокореляційної функції, старший індекс яких приймається за порядок авторегресійної моделі.

При цьому враховують, що вибіркові оцінки автокореляційної функції наближено [1] характеризуються зміщенням і дисперсією (середньоквадратичне відхилення ). Якщо значення не потрапляють до довірчого інтервалу , тобто , вони вважаються значущими. У протилежному випадку, тобто коли

, (71)

приймаємо гіпотезу .

На рис. 4.3 штриховими лініями позначено область значущих значень оцінок :

Рис.6. Інтервали значущих оцінок автокореляційної функції

Після орієнтовного визначення структури регресійної моделі для її подальшого уточнення необхідно виконати оцінювання параметрів авторегресії. Зокрема для моделі виду (65) це може бути реалізоване методом найменших квадратів (МНК). Приймаючи до уваги структуру моделі (65), сформуємо із масиву нев'язок вектор зашумлених значень залежної змінної Z та матрицю плану Х:

; . (72)

Обробка цих даних за МНК дає вектор МНК-оцінок параметрів авторегресійної моделі.

Однак як вже зазначалося, введена орієнтовна оцінка р порядку авторегресії може значно перевищувати фактичне значення порядку, а структура моделі - містити надлишкові регресори. Тому необхідно обов'язково провести перевірку та уточнення структури моделі шляхом виявлення за критерієм Стьюдента незначущих параметрів та виключення із структури відповідних їм регресорів (див. п.5.3).

Кількість регресорів, що звичайно має бути вилучена з вихідної структури, може значно перевищувати число регресорів, які залишилися. Порядок вилучення припускає високий ступінь суб'єктивізму при прийнятті проміжних рішень. Через це результатом спрощення може бути кілька достатньо суперечливих кінцевих варіантів структури. Тому для визначення структури авторегресії рекомендується використовувати суттєво різні підходи [1,2,12]. Для отримання кращої структури моделі результати, отримані за різними підходами, треба співставити та проаналізувати. Викладений вище метод виключення регресорів із початкової структури моделі доцільно комбінувати з методами включення регресорів. Це вже розглянутий у п.3 метод крокової регресії, а також евристичний підхід, в основі якого включення до початкової структури моделі лагових змінних, лаг яких співпадає з індексами вибіркових коефіцієнтів ,..., значення яких виходять за межі обрахованого вище 95% довірчого інтервалу для оцінок функції автокореляції (при чому для індексів виділених вибіркових коефіцієнтів має виконуватися умова: ). Методи включення дозволяють ефективно контролювати процес ускладнення початкової моделі за допомогою критерію Фішера та коефіцієнта детермінації R2, комбінуючи цей контроль із перевіркою значимості коефіцієнтів моделі за Стьюдентом. Цей процес триває доки не буде отримана модель, адекватна вихідним даним, яка і буде використана для коригування оцінок тренду. В загальному випадку, деталізуючи наведене вище співвідношення (66), що описує спосіб введення корекції тренду, з урахуванням моделі (65) отримуємо:

 (73)

Об'єктивним показником ефективності використання побудованої моделі авторегресії має бути некорельованість значень нев'язки , розрахованої за скорегованими оцінками тренду:

, . (74)

Тобто вибіркові коефіцієнти автокореляційної функції r1, r2, ..., розраховані за рядом , мають перебувати в межах введеного вище 95% довірчого інтервалу.

9. Прогноз детермінованої складової часового ряду. Метод експоненціального згладжування.

9.1 Локальні моделі часового ряду

Відповідно до змісту адитивної моделі часового ряду (40), детермінована складова цієї моделі в загальному випадку відрізняється від тренду , включаючи в себе інші компоненти, наявність яких обумовлює необхідність застосування для опису складової певних локальних моделей. Особливість цих моделей полягає в тому, що в ряді потреби вони можуть змінювати свою структуру та параметри, адаптуючись до появи у складі додаткових компонентів (зокрема циклічного компоненту ).

Припустимо, що інформативна складова , моделі (40) часового ряду отримана в результаті регулярної дискретизації процесу , який в околі точки t може бути розвинений у ряд Тейлора

, (75)

причому для цілком достатньої точності наближення можна обмежитися першими k+1 членами ряду, тобто поліном

. (76)

Цей поліном , коефіцієнти якого визначені із врахуванням особливостей поведінки функції в околі точки t, являє собою локальну модель часового ряду. Точкою локалізації цієї моделі є рівень (або t -й момент часу). В принципі точкою локалізації може бути будь-яке значення часового ряду , , що певним чином наближає локальну модель до моделі зваженого ковзкого середнього. Однак модель ковзкого середнього будується лише за даними, що входять до складу КІ, тоді як локальна модель - за даними фрагменту часового ряду , причому якщо , то цей фрагмент - весь часовий ряд.

Одним з основних питань побудови локальної моделі є визначення вектора її коефіцієнтів . Ця задача розв'язується шляхом застосування зваженого МНК. Зокрема, у випадку коефіцієнти поліному визначаються з умови мінімуму зваженої суми квадратів нев'язок виду:

, (77)

де , . В результаті отримуємо локальну модель часового ряду, яка добре адаптується до його особливостей в околі точки , виказуючи високі апроксимативні якості в цьому околі, однак поступово втрачає точність апроксимації з віддаленням від . Цими властивостями модель завдячує вектору своїх коефіцієнтів , визначеному за критерієм (77), точніше системі ваг , :

(78)

значення яких спадають за експоненціальним законом так, що квадрати нев'язок в околі мають найвищу вагу (чим й забезпечується найкраща апроксимація даних цієї області ряду), яка потім монотонно зменшується, наближаючись до нуля для початкових індексів елементів часового ряду.

Якщо до локальної моделі (76), побудованої для , ввести значення змінної , отримаємо прогноз часового ряду на інтервал упередження L. Вибір значення L залежить від ширини околу, в якому локальна модель задовільно апроксимує часовий ряд. Окіл ширшає із збільшенням порядку k та зменшенням значення і навпаки, вужчає із зростанням та зменшенням k. Звичайно , а , однак іноді ліва межа для параметру може спадати навіть до 0,3.

У разі збільшення обсягу часового ряду (звичайно n зростає в результаті появи нових даних, що поповнюють часовий ряд справа), локальна модель має постійно адаптуватися до нових даних, кожного ряду змінюючи значення вектора коефіцієнтів відповідно до постійно обновлюваної статистики в лівій частині критерію (77).

Загалом наведена вище процедура застосування локальної моделі для прогнозу часового ряду видається достатньо копіткою, особливо в разі покрокового оновлення ряду за рахунок поелементного приєднання до нього справа нових даних. В цьому випадку виникає необхідність в багаторазовому застосуванні зваженого МНК для перерахунків значень вектора на покроково зростаючій вибірці даних. Тому видається доцільним застосування рекурентних методів обробки часового ряду, які дозволяютьє суттєво спростити обчислення коефіцієнтів локальних моделей , . Найбільш поширеним серед них можна вважати метод експоненціального згладжування - модифікацію дисконтованого методу зважених найменших квадратів, пристосовану для рекурентного прогнозу елементів часового ряду.

9.2 Експоненціальне згладжування

Ідея методу експоненціального згладжування [1,3] полягає в згладжуванні часового ряду ковзкою середньою з експоненціальними вагами (78). Така середня більше характеризує значення процесу на кінці інтервалу згладжування, ніж на початку. Справді, здебільшого на прогнозні значення суттєво впливають останні рівні часового ряду, тому їм надається більша вага порівняно з початковими, які, однак, не виключають зовсім з аналізу, тому що вони також несуть певну інформацію про досліджуваний процес. Назва методу випливає з того, що дані згладжуються за допомогою зваженої середньої, у якій ваги змінюються згідно з експоненціальним законом.

Ряд спостережень можна згладжувати будь-яку кількість разів, наприклад, ряд, згладжений один раз, можна згладити ще один раз ковзкою середньою з експоненціальними вагами (так зване подвійне експоненціальне згладжування), отриманий ряд знову згладжуємо ковзкою середньою з експоненціальними вагами (потрійне експоненціальне згладжування) тощо.

Метою багаторазового експоненціального згладжування є рекурентний розрахунок оцінок коефіцієнтів рівняння виду (75), який базується на обчисленні експоненціальних середніх (ЕС) різних порядків.

ЕС першого порядку розраховують за вихідною вибіркою :

, (79)

де - постійна згладжування (звичайно обирається у діапазоні , але іноді може зростати до ), - фактор затухання.

Використання рекурентної формули (79) потребує завдання початкового значення при (див. нижче підрозділ "Вибір початкових умов").

Покажемо, як залежать ЕС першого порядку (формула (79)), обрахована на момент часу , від раніш отриманих згладжених величин:

Як видно, є лінійною комбінацією всіх попередніх рівнів, ваги яких зменшуються в геометричній прогресії. Вага рівня, віддаленого на відліків від поточного моменту , дорівнює . Наприклад, якщо постійна згладжування , то для моменту часу:

рівень буде мати вагу ,

-"- ,

-"- ,

-"- тощо.

ЕС другого порядку розраховують за вибіркою ЕС першого порядку за формулою:

, (80)

аналогічно ЕС третього порядку - за вибіркою ЕС другого порядку:

. (81)

Для практичного обчислення ЕС можна скористатися вбудованою можливістю Excel (Сервіс Аналіз даних Експоненціальне згладжування).

9.3 Розрахунок оцінок параметрів і прогнозних значень

Обрахувавши експоненціальні середні різних порядків, далі обчислюються оцінки параметрів у виразі (75), та прогнозні значення рівнів ряду. Співвідношення, за якими виконуються обчислення, залежать від того, поліном якого ступеня використовується у локальній моделі часового ряду. Зазначимо, що ступінь апроксимуючого поліному (75) звичайно не перевищує 2.

У довільний момент часу прогноз на відліків уперед обчислюється (як це випливає з виразу (76)) за формулою:

. (82)

Для поліному степеня , тобто для простої середньої, ця формула трансформується у вираз:

. (83)

Прогноз часового ряду на точок вперед лінійним поліномом () розраховується за формулою:

, (84)

, . (85)

Квадратична поліноміальна модель () для прогнозування має вигляд:

, (86)

. (87)

9.4 Вибір початкових умов

Для розрахунку ЕС необхідно задати початкову умову при . Єдиного підходу до завдання початкових наближень немає, часто їх задають відповідно до змісту даних чи умов дослідження. Наприклад, як початкове наближення обирають , особливо, коли аналізується ряд спостережень з флуктуаціями (біржовий курс тощо). Інколи обирають середнє значення, особливо коли даний ряд коливається біля свого середнього (продаж товару, попит на який стабільний у часі тощо): . Крім того, як початкові наближення можна обирати середні декількох початкових значень.

Якщо використовувати початкові умови першого варіанта, то при згладжене значення попереднього етапу (тобто початкове наближення) приймається рівним першому значенню вихідного ряду , тоді:

. (88)

Для наступних моментів часу використовують формулу (79), наприклад:

, ,

, ,

, .

9.5 Вибір оптимальних параметра згладжування та степеня полінома методом ретропрогнозу

Другою проблемою при прогнозуванні методами експоненціального згладжування є вибір параметра згладжування . Слід врахувати, що значення параметра суттєво впливає на результат прогнозу. Якщо значення близько до одиниці, то при прогнозуванні враховують в основному вплив останніх рівнів ряду; якщо близьке до нуля, то ваги, за якими зважуються рівні ряду, спадають повільно, що дає можливість врахувати практично всі попередні значення ряду.

Виконуючи прогноз, слід уточнити степінь поліноміальної моделі та параметр згладжування . Для кожних обраних значень та будуються моделі прогнозу на відліки та обчислюють прогноз на останніх рівнях ряду, тобто виконують ретропрогноз (якщо обсяг вихідних даних дозволяє, бажано брати ). Для цих рівнів розраховують відхилення значень ретропрогнозу за локальною моделлю від фактичних даних часового ряду та дисперсії цих відхилень за формулою:

. (89)

Уточнення виконують шляхом зіставлення дисперсій відповідних прогнозів для кожного параметра згладжування , степеня поліному при прогнозі на точки.

Двійка значень , для яких буде отримано мінімум показника (окремо для кожного інтервалу прогнозу ), вважаються оптимальними й використовуються для подальших розрахунків прогнозу за локальною моделлю .

10. Особливості згладжування часових рядів в присутності аномальних даних. Усунення аномальних даних

10.1 Виявлення аномальних даних в часових рядах

При попередньому аналізі вихідних даних іноді можна спостерігати окремі значення, які суттєво відрізняються від значень, що знаходяться поруч. Це так звані аномальні дані, що виникли, наприклад, через випадкові збої у вимірювальній або обчислювальній техніці. Аномальні дані (АД) можуть зустрічатися у будь-яких експериментальних даних: як у вибіркових сукупностях, так і у часових рядах, причому у кожному з цих випадків не вилучені АД стають причиною появи серйозних помилок у результатах , отриманих з використанням цих даних.

Для розробки методів виявлення АД в часовому ряді (послідовності) , принциповим є питання сталості, незмінності в часі імовірнісних характеристик процесу . Якщо це має місце, процес називається стаціонарним у вузькому сенсі. Для такого процесу розподіл ймовірностей його значень лишається однаковим і постійним для будь-якого моменту часу . За цих умов до виявлення АД серед елементів часового ряду , можна залучати методи виявлення АД у вибіркових сукупностях, моделлю яких є множина значень, отриманих у серії випробувань випадкової величини . Однак зазначені методи не враховують існуючих між елементами часового ряду кореляційних зв'язків, наявність яких дозволяє зробити пошук та виокремлення АД більш успішним.

При наявності відомостей про автокореляційну функцію процесу доцільно спершу від вихідного часового ряду перейти до різницевого ряду ,, де

 . (90)

За відсутності АД математичне сподівання такого ряду дорівнює нулю, а дисперсія може стати меншою за :

, (91)

тобто за сильної автокореляції для виявлення АД більш сприятливі умови будуть при обробці різницевого ряду. Слід також враховувати, що поряд з можливим зменшенням дисперсії відносно при наявності АД після переходу до різницевого ряду виникає своєрідне "підсилення" аномальності за рахунок розмноження АД. Наприклад, якщо елемент є аномальним, аномальність (різного знаку) будуть мати вже дві різниці: та .

Якщо ряд нестаціонарний, розв'язок задачі пошуку та виявлення АД суттєво ускладнюється і стає залежним від характеру нестаціонарності.

Коли ряд нестаціонарний за математичним сподіванням і структуру його елемента можна представити у вигляді:

, (92)

є можливість шляхом згладжування вихідного ряду оцінити трендову складову, обчислити нев'язки і далі, аналізуючи послідовність нев'язок, виявити АД. Звичайно уживані методи виділення тренду шляхом поліноміальної апроксимації вихідного ряду або лінійних згладжуючих фільтрів при наявності у вихідній послідовності АД приводять до надто суттєвих зміщень в отриманих оцінках тренду. Через це перспективи подальшої роботи з цими оцінками вельми суперечливі, тому для оцінювання тренду слід застосовувати робастні алгоритми, серед яких найбільш поширеним є медіанне згладжування.

10.2 Метод медіанного згладжування

Застосовування методу медіанного згладжування для виявлення та усунення АД реалізується шляхом аналізу даних в межах ковзкого інтервалу, що покроково переміщується вздовж часового ряду в заданому напрямі, звичайно за наростанням індексів членів ряду. Довільному положенню ковзкого інтервалу відповідає фрагмент вихідного часового ряду з елементів: ,…,,…,. Згладжена оцінка визначається для середнього елемента інтервалу шляхом побудови варіаційного ряду з його елементів ,…,,… та виділення його медіани:

. (93)

Після цього зі складу ковзкого інтервалу виключається крайній лівий елемент та дописується справа черговий елемент часового ряду. При цьому середина інтервалу зміщується на один елемент вправо. Для нового положення ковзкого інтервалу повторюється процедура обчислення згладженої оцінки середнього елемента. У такий спосіб ковзкий інтервал проходить повздовж усього часового ряду, починаючи з його лівого краю.

Після згладжування вихідного ряду можна оцінити трендову складову часового ряду і далі, аналізуючи послідовність нев'язок, виявити АД.

Якщо у вихідній вибірці відсутні АД кратності (кратність - кількість послідовних елементів ряду, що являють собою аномальності), то в якості згладжених оцінок елементів часового ряду будуть виступати тільки кондиційні елементи ряду, які не належать до аномалій. Таким чином, якщо у вихідному часовому ряді присутні АД, кратність яких не перевищує , в згладженому ряді їх не буде, проте в часовому ряді, утвореному нев'язками, будуть присутні відповідні АД з вихідного ряду.

З точки зору забезпечення надійного робастного згладжування, стійкого до АД високої кратності , слід використовувати ковзкий інтервал з великим . З іншого боку, із ростом , за рахунок так званої динамічної похибки медіанного фільтра, починають зростати амплітуди нев'язки для кондиційних даних, що стає завадою для надійного виділення АД на фоні збільшених (зрослих) значень нев'язок.

Одним з рішень в цьому випадку є подвійне медіанне згладжування, в якому послідовність згладжених даних , отриманих на першому етапі обробки вихідного часового ряду (з достатньо великим ), ще раз «підгладжується» медіанним фільтром з . Прикладом такої обробки є так звана процедура „Тьюки 53”, де , . Отримавши після другого згладжування ряд , фінальну згладжену оцінку рекомендується розраховувати за формулою:

. (94)

Більш складні випадки нестаціонарності часового ряду вимагають застосування більш складних методів обробки даних.

Література

1. Френкель А.А. Прогнозирование производительности труда: методы и модели. - М.: Экономика, 1989. - 214 с.

...