Регресійний аналіз. Множинна регресія. Лінійна,

поліномінальна та логістична регресії. Кластерний аналіз.

Стандартизація даних. Ієрархічні та дивизимні методи.

Міри відстані та схожості, правила об’єднання. Аналіз

отриманих результатів.

матеріали використано з доробку

Кафедра теорії ймовірностей, статистики

та актуарної математики, Київського

національного університету

імені Тараса Шевченка

Лінійна регресія

Створимо новий файл, в якому змінну VAR1 заповнимо послідовно

значеннями від 0 до 10, змінну VAR2 - випадковими значеннями від 0 до 1,

а змінну VAR3 задамо, як суму VAR1+ VAR2.

Виконаємо послідовність команд: Statistics Æ Basic Statistics/Tables Æ

Descriptive Statistics Æ Prob.&Scatterplots (див. рис. 1).

Рис. 1

Як Variables виберемо VAR1-VAR3, натиснемо кнопку 2D Scatterplot.

У першому списку змінних вкажемо VAR1, в другому

- VAR3 і

натиснемо кнопку ОK (див. рис. 2).

На графіку, що з’явився, зображено пряму лінійної регресійної моделі

для VAR3 через VAR1, а у верхній частині вікна бачимо рівняння лінійної

регресії (див. рис. 3).

Якщо у змінній VAR3 замінити одне із значень, наприклад на 70, і

побудувати графік знову, то побачимо, що рівняння регресії буде

враховувати дане значення і з графіка буде очевидно, що 70 є викидом

(див. рис. 4).

Натиснемо піктограму Brushing. У вікні, що з’явиться зробимо

активними Exclude та Box, виділимо прямокутником значення (див. рис.

5) і натиснемо кнопку Apply. Виділене значення зникне з графіка і

регресійна пряма змінить своє положення.

Для того, щоб значення викиду не виводилось у наступних графіках та

не враховувалось при обчисленні регресійної формули, у змінній VAR4

заповнимо всі значення одиницями, а те значення, що стоїть напроти

викиду - нулем (див. рис. 6).

Рис. 4

У вікні Descriptive Statistics натиснемо кнопку Weight (див. рис. 1).

Оберемо змінну VAR4, перемкнемо Status на On та натиснемо OK (див.

рис. 7). Тепер при аналізі викиди враховуватися не будуть.

Рис. 7

Нехай маємо таку таблицю з даними:

100

200

300

400

500

600

700

де Y - врожайність, X - добрива, Z - опади.

Потрібно знайти формулу багатофакторної лінійної регресії для Y:

Y=B

X+B Z,

де B

- невідомі коефіцієнти.

Виконаємо послідовність команд: Statistics Æ Multiple Regression, як

змінні оберемо Y - залежна, X і Z - незалежні (див. рис. 8). Натиснемо

OK. Отримуємо результат, який зображено на рис.

В закладці Quick натиснемо кнопку Summary: Regression results.У вікні,

що з’явилось (див. рис.

10), бачимо оцінки параметрів та допоміжну

статистику. Обидві змінні значимі (виділені червоним). У третьому стовпці

таблиці вказані оцінки для коефіцієнтів B . Отже,

Y=2095+0.038*X+0.833*Z.

Якщо наші змінні попередньо стандартизувати, то у результаті такого

регресійного аналізу отримали б оцінки коефіцієнтів, які записані в

першому стовпчику таблиці (зрозуміло, що В₀=0). Коефіцієнти з першого

стовпця показують внесок у регресійну модель змінних X та Z.

Рис. 10

Для того, щоб обчислити передбачуване значення для Y для заданих X

та Z і побудувати

95% проміжок надійності, перейдемо у закладку

Residuals/assumptions/prediction і натиснемо Predict dependent variable (див.

рис. 11). У відповідні віконця вводимо значення змінних X та Z (див. рис.

12) і натискаємо OK.

Якщо незалежні змінні набувають одного і того ж значення його можна

ввести у віконці Common Value і натиснути Apply.

У вікні результатів аналізу (див. рис.

13) бачимо передбачуване

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

Рис. 13

Якщо ж потрібно подивитись як розподілені залишки, то натиснемо

кнопку Perform residual analysis. У вікні, що з’явилось (див. рис.

14),

зібрані різні методи для аналізу залишків регресійної моделі.

Рис. 14

Наприклад, натиснувши Normal plot of residuals отримаємо Q-Q графік,

на якому видно наскільки залишки узгоджуються з нормальним законом

розподілу (див. рис. 15).

Багатофакторна регресія

Відкриємо файл Job_prof.sta

(див. рис.

1). У файлі вказано бали,

отримані претендентами під час тестування при прийнятті на посаду (у

перших чотирьох стовпцях), та оцінка професійної здатності претендентів

(у п’ятому стовпці) після закінчення випробувального терміна. Нам

потрібно знайти лінійну багатофакторну регресійну модель залежності

оцінки професійної здатності від оцінок за тести. Завантажимо модуль

Multiple Regression: Statistics Æ Multiple Regression.

Рис. 1

Натиснувши Variables, обираємо Job-Prof, як залежну змінну, а як

незалежні змінні вибираємо перші чотири змінні (див. рис. 2). Двічі

натискаємо ОK. Результати регресійного аналізу зображені на рисунку 3.

Всі змінні, окрім другої, є значимими (виділені червоним).

У закладці Quick натиснемо Summary: Regression result (див. рис. 3). У

вікні, що з’явилось (див. рис. 4), бачимо результати аналізу: у третьому

стовпці - коефіцієнти багатофакторної лінійної регресійної моделі, а в

першому стовпці

- коефіцієнти цієї ж регресійної моделі для

стандартизованих змінних.

Проаналізувавши результати прийдемо до висновку, що Test 2 досить

мало впливає на оцінку професійної здатності претендентів: відповідний

коефіцієнт у першому стовбці становить 0,043. Тому можливо є доречним

взагалі вилучити Test 2 з регресійної моделі.

Рис. 2

Проведемо спочатку аналіз залишків побудованої моделі. У закладці

Residuals/assumptions/prediction, натиснемо кнопку Perform residual

analysis. У вікні, що з’явилось, обираємо Residuals Æ Casewise plot of

residuals (див. рис. 5).

Рис. 5

Так ми перевіримо, чи не виходять залишки за межі 3δ. Отримаємо

таблицю (див. рис. 6), в якій знаком «*» вказано де знаходиться залишок

у інтервалі

(-3δ;

3δ). У правій частині таблиці міститься додаткова

інформація про залишки.

Бачимо, що залишки лежать у проміжку (-2δ; 2δ) і їхні середнє і медіана

дорівнюють нулеві.

Потім оберемо закладку Outliers (див. рис. 7). Натиснемо Casewise plot

outliers та упевнимося, що викидів немає

- з’явиться відповідне

повідомлення (див. рис. 8).

На основі цих результатів можна вважати, що багатофакторна

регресійна модель достатньо добре описує наші дані.

Цікаво дослідити, що відбудеться, якщо вилучити змінну Test 2, яка

досить мало впливає на оцінку професійної здатності претендентів, з

регресійної моделі. Розглянемо, як автоматизовано процес знаходження

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

У вікні Multiple Regression у закладці Advanced відмітимо Advan-ced

options (stepwise or ridge regression) та натиснемо OK (див. рис. 9).

Рис. 9

У закладці Stepwise вкажемо метод Forward stepwise - змінні будуть

введені у регресійну модель по одній. У Display results вкажемо

покроковий вивід результатів

- At each step

(див. рис.

10). Тобто ми

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

кожного кроку. У полі F to enter вказуємо значення 1, а у полі F to remove

вказуємо 0.01. Ці два числа визначають верхню та нижню межі проміжку

для значимості внеску у регресійну модель змінних. Якщо значимість

змінної потрапляє в цей проміжок, то включаємо її до регресійної моделі,

інакше - відкидаємо.

Рис. 10

Оскільки маємо чотири незалежні змінні, то кількість кроків множинної

регресії Number of steps досить вказати рівною чотирьом (після кожного

кроку до моделі може включатись не більше однієї змінної). Натискаємо

кнопку OK.

У вікні, що з’явилось, бачимо, що в моделі ще немає жодної змінної

(див. рис. 11). Натискаємо Next. З’явився перший коефіцієнт та перша,

вибрана до рівняння регресії змінна (див. рис. 12). Натискаємо кнопку

Next, допоки на цій кнопці не з’явиться напис ОК. Це означатиме, що

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

значимість внеску яких знаходиться у вказаних межах увійшли до рівняння

регресії. Бачимо, що змінна Test 2 не увійшла у нову багатофакторну

регресійну модель (див. рис. 13).

Далі можемо подивитися результати. Для цього натиснемо на закладці

Advanced, кнопку Stepwise Regression Summary (див. рис. 14). У вікні, що

з’явилось (див. рис. 15), бачимо, статистику внеску обраних змінних і

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

Advanced, кнопку Summary: Regression results (див. рис. 14), то у вікні, що

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

до регресійної моделі.

Використовуючи отриману інформацію

(рис.

15 і

16) можемо

порівняти стару і нову регресійну моделі. Бачимо, що вони майже не

відрізняються. Отже, внесок змінної Test 2 дійсно був незначним.

Рис. 15

Рис. 16

Якщо є мультиколінеарність змінних (наприклад, в кількох стовпцях

містяться дані про ціну одного і того ж товару в різних валютах), то

очевидно, що для регресійної моделі слід взяти не всі із таких змінних.

Якщо для обробки таблиці мультиколінеарних даних так як і раніше

скористатися стандартним алгоритмом, то з’явиться повідомлення про

помилку (див. рис. 17).

Рис. 17

У цьому разі для регресійного аналізу потрібно вибирати метод

Forward stepwise (який дозволяє вводити у модель по одній змінній), а не

Standard.

Нелінійна регресія

Досить часто із діаграми розсіювання, або попереднього досвіду роботи

з подібними даними зрозуміло, що регресійна модель не є лінійною. Отже,

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

входять не просто змінні, а деякі функції від них (наприклад, степінь,

експонента, логарифм). Покажемо, як діяти у цьому випадку.

Створимо новий файл. Першу змінну заповнимо числами від 1 до 10,

другу - випадковими значеннями від 0 до 1, третя змінна

v3=2*v1*v1+v2+20.

Потрібно отримати рівняння регресії для VAR3 через VAR1. Зрозуміло,

що доцільно ввести у регресійну модель квадрат VAR1.

Виконаємо Statistics Æ Advanced Linear/Nonlinear Models Æ Fixed

nonlinear Regression (див. рис. 1).

Рис. 1

У вікні, що з’явиться

(див. рис.

2), оберемо зміні VAR1 та VAR3.

Натискаємо OK. Оскільки ми знаємо, що наша залежність квадратична, то

виберемо X**2 (див. рис. 3). Натискаємо OK.

З’явиться вікно Model Definition. Перейдемо на закладку Quick і

вкажемо залежну змінну - VAR3, а незалежними оберемо змінні VAR1 та

VAR1**2 (див. рис. 4).

Рис. 4

Будуємо регресійну модель покроково з відкиданнями незначимих

змінних. У закладці Stepwise вказуємо значення 1.0 для F to enter, а у полі F

to remove вказуємо 0.01. Як і очікували, у регресійній моделі одержали

єдину значиму змінну VAR1**2 (див. рис. 5).

Натиснувши Summary: Regression results отримуємо таблицю

результатів регресійного аналізу

(див. рис.

6). Оцінки коефіцієнтів

регресії близькі до значень у заданій теоретичній моделі.

xλ,

≠

ln( x), λ

≠

Використаємо для цього допоміжну програму. Вона обчислить оцінку

максимальної вірогідності параметра λ. Для цього виконаємо

послідовність команд Help Æ Open Examples Æ Macros Æ Analysis

Examples і викличемо програму BoxCox. У вікні, що з’явилось (див. рис.

7), бачимо текст програми на мові STATISTICA Visual Basic.

Рис. 7

Натиснемо на стрілку в панелі макросів, щоб почати виконання

програми. Вказуємо наші змінні і у вікні, яке потім з’являється (див. рис.

8), задаємо межі для параметра λ . Натискаємо OK. Вибравши один із

запропонованих методів

(див. рис.

9) Manual/visual search

(graph of

lambda vs. SSE) або Iterative optimization(golden search) отримаємо оцінку

для λ, яку потім можна використати при побудові регресійної моделі.

Рис. 8

Рис. 9

За допомогою програми Boxtidwell.stb, яка знаходиться в тій же

директорії, що і BoxCox,

можна обчислити оцінки максимальної

вірогідності для параметрів перетворення незалежних змінних.

Розглянемо деякі інші нелінійні регресійні моделі.

Відкриємо файл program.sta

(див. рис. 10). У першій колонці файлу

вказано витрачений працівниками час на вивчення певного програмного

продукту, а у другій колонці вказано результати тестового завдання, на

виконання якого було відведено фіксований час. Потрібно побудувати

залежність успішності виконання від часу вивчення. Шукатимемо цю

залежність у вигляді логістичної регресійної моделі:

exp(b

1+exp(b

де

- успішність, t - час,

b₁

- невідомі параметри.

Завантажимо модуль з нелінійного оцінювання: Statistics Æ

Advanced/Linear Nonlinear Models Æ Nonlinear Estimation (див. рис. 11).

Обираємо Quick Logit regression (див. рис. 12) і натискаємо OK.

Рис. 12

Вкажемо як залежну змінну Success, а як незалежну - Experience (див.

рис. 13). Натискаємо двічі OK. Отримані результати відобразяться у вікні,

яке зображено на рисунку

Рис. 13

Натиснувши Fitted 2D function & observed values отримаємо графік

логістичної регресійної функції (див. рис. 15), яка найкраще описує наші

дані. Формула, яка визначає функцію, написана над графіком. На цьому ж

графіку зображено наші дані. Їхні координати по вертикальній осі є або 1,

якщо тестове завдання виконано успішно, або 0, якщо ні. Використовуючи

отриману логістичну модель можна передбачити наскільки успішно

Натискаючи різні кнопки у закладці Residuals (див. рис. 16), можемо

отримати гістограму залишків, значення самих залишків, значення оцінок,

які

«пророкують» згідно із затраченим часом, вивчити відповідність

розподілу залишків до гауссового розподілу. Наприклад і Q-Q графік і

гістограма

(див. рис.

17 і

18) свідчать про нормальний розподіл

залишків.

Рис. 16

Рис. 17

Рис. 1

Відкриємо файл Learning.sta

(див. рис.

19).

Цей файл містить

інформацію, про два підприємства А і В, які переходять на випуск нової

продукції. Одна змінна - час, інша ефективність - відношення прибутку на

один виріб нової моделі до старої. З часом на підприємствах все краще

освоюють технологію пов’язану з випуском нової продукції, тому

ефективність зростає. Ми хочемо порівняти два підприємства між собою і

визначити яке підприємство має більший потенціал для впровадження

і моделі для даних

нових технологій. Для цього ми побудуємо регресійн

кожного підприємства і порівняємо ці моделі.

Виконаємо Statistics Æ Advanced/Linear Nonlinear Models Æ Nonlinear

Estimation Æ User specified regression, custom loss function

(див. рис. 12).

Як Function to be estimated & loss function задамо функцію

v3=B

*V1+B

*exp(B

*V2),

а як міру відхилення (OBS-PRED)**2 (див. рис. 20). Натискаємо ОК.

Оскільки текстовому значенню першої змінної Plant_A відповідає 1, а

текстовому значенню Plant_B відповідає 0, то додатний знак коефіцієнта В

У вікні Results, що з’явилось (див. рис. 22), бачимо результати аналізу.

Значення Proportion of variance accounted for показує наскільки пояснено

розкид даних.

Побачити оцінки для параметри моделі можемо, натиснувши кнопку

Summary (див. рис.

23). Оскільки В <

0, то друге підприємство має

більший потенціал для впровадження нових технологій. Це очевидно, якщо

подивитись на малюнок 19 з даними, які ми використовували для нашого

дослідження.

Аналіз залишків, як і раніше, можна здійснити використовуючи меню з

закладки Residuals (див. рис. 22). У результаті приходимо до висновку,

що модель добре описує наявні дані.

Рис. 22

Рис. 23

Часові ряди (виділення

періодичних складових)

Створимо нову таблицю даних, у якій 100 спостережень (Cases). Першу

змінну VAR1 задамо формулою:

= 0.9*sin(2*Pi*0.4*V0) + 0.4*sin(2*Pi*0.2*V0) + 4

де V0 - це змінна, що послідовно набуває значень від 1 до 100 (тобто

номер спостереження по порядку). Зауважимо, що змінна VAR1 є сумою

двох періодичних функцій. Числа 0.4 та 0.2 - це частоти коливань, а числа

2.5=1/0.4, 5=1/0.2, відповідно, періоди коливань.

0.9 та 0.4 - амплітуди

першого та другого коливань.

Проаналізуємо цей часовий ряд. Виконаємо Statistics Æ Advanced

Linear/ Nonlinear Models Æ Time Series/Forecasting (див. рис. 1).

Рис. 1

У вікні, що з’явиться (див. рис. 2), виберемо у Variables змінну VAR1 і

натиснемо кнопку Spectral (Fourier) Analysis. Для того, щоб побачити

графік нашого ряду, у вікні що з’явиться перейдемо у закладку Review

series і натиснемо першу кнопку Plot біля кнопки Review highlighted

variables (див. рис. 3).

Натиснемо кнопку Periodogram, у вікні що з’явилось (див. рис. 6).

З’явиться графік

(див. рис.

7), піки якого відповідають частотам

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

випадку це 0.4 та 0.2.

Кнопка Spectral density дає змогу побудувати згладжену періодограму за

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

Advanced. Оберемо, наприклад, Daniel (див. рис. 8). Отримуємо графік

(див. рис. 9), який відрізняються від періодограми з малюка

Змінюючи внизу вікна опції з Frequency на Period і натискаючи кнопку

Periodogram отримуємо графік (див. рис. 10) з піками, що відповідають

періодам

2.5 та

5. Таким чином можемо визначити найважливіші

періодичні компоненти, що присутні у змінній VAR1.

Рис. 6

Натиснувши кнопку Summary (див. рис. 6, 8) отримаємо таблицю

(див. рис. 11), у якій перші два стовпці - це частоти і періоди, які їм

відповідають. У наступник стовпцях показано з яким коефіцієнтом

відповідні синуси чи косинуси входить у ряд Фур’є розкладу змінної VAR1.

Червоним кольором виділені великі значення, що відповідають частотам

0.2 та 0.4; решта коефіцієнтів близькі до нуля.

Рис. 11

Відкриємо файл Sunspot.sta. У цьому файлі міститься щорічна

інформація про кількість плям на сонці. Будемо вивчати періодичність

сонячної активності.

Виконаємо Statistics Æ Advanced Linear/ Nonlinear Models Æ Time

Series/Forecasting Æ Spectral

(Fourier) Analysis. Побудуємо, як і в

попередньому прикладі, графік часового ряду

(див. рис. 12). Бачимо, що

у даних спостерігається певна періодичність.

Рис. 12

Перейдемо у закладку Quick і оберемо Single Fourier Analysis.

Відмітимо опцію Frequency і натиснемо кнопку Spectral density. На графіку,

що з’явився (див. рис.

13), спостерігаємо дві періодичні складові, які

домінують всі інші.

Натиснувши кнопку Summary отримаємо таблицю з коефіцієнтами

синусів та косинусів для кожної з частот. Проте з даної таблиці важко

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

одразу побачити яким частотам відповідають найбільші за модулем

коефіцієнти, перейдемо у закладку Advanced і натиснемо кнопку N largest

values. В отриману таблицю

(див. рис.

14) виведено

10 значень з

найбільшими амплітудами. Бачимо, що основний період сонячної

активності становить трохи більше 11 років. Є ще інша вагома періодична

компонента з періодом близько 80 років.

Якщо у закладці Quick відмітити Period і натиснути Periodogram, то на

графіку, що з’явився

(див. рис.

15), спостерігаємо ці дві періодичні

складові, які домінують всі інші.

Рис. 13

Часові ряди

(кореляційний аналіз)

Відкриємо файл Stocks.sta. У ньому є дані про коливання щоденних

курсів акцій двох типів. Виконаємо Statistics Æ Advanced Linear/ Nonlinear

Models Æ Time Series/Forecasting. Як Variables виберемо першу змінну

Stock1 (див. рис. 1).

Рис. 1

Для аналізу часового ряду натиснемо кнопку ARIMA & autocorrelation

functions. Перейдемо на закладку Review series і натиснемо кнопку Plot біля

Review highlighted variables (див. рис. 2). На графіку ряду

(див. рис.

3) бачимо, що навколо тренда присутні високочастотні коливання. Для

подальшого аналізу ряд спочатку потрібно згладити.

Перейдемо на закладку Advanced і натиснемо кнопку Other

transformations & plots

(див. рис. 4). У вікні Transformations of Variables

виберемо у закладці Smoothing опцію N-pts moving average N=5 (тобто

заміняємо кожне значення у часовому ряді на середнє арифметичне п’яти

Рис. 6

Для того, щоб порівняти вихідний ряд та щойно отриманий, перейдемо

у закладку Review & plot і натиснемо кнопку Plot біля Review multiple

variables

(див. рис. 7). У вікні Select Variables for the Spreadsheet/Plot,

що з’явилось

(див. рис.

8), оберемо імена змінних до і після

перетворення. Натискнемо OK. На отриманому графіку

(див. рис.

різними кольорами зображено початковий часовий ряд і ряд після

згладжування.

Якщо ми хочемо зберегти модифікований ряд, то можемо натиснути

Save variables

(див. рис.

5). З’явиться таблиця зі значеннями часового

ряду до і після перетворення.

Рис. 9

Повернімося у вікно Time Series Analysis. Як Variables виберемо Stock2.

Натиснемо ARIMA & autocorrelation functions. Перейдемо на закладку

Review series і натиснемо кнопку Plot біля Review highlighted variables (див.

рис. 2). Отримуємо графік ряду (див. рис. 10).

Виділимо з ряду лінійний тренд. Для цього перейдемо на закладку

Advanced і натиснемо кнопку Other transformations & plots

(див. рис. 4).

У панелі що з’явиться, перейдемо на закладку x=f(x), виберемо Trend

subtract (x=x-(a+b*t)) (див. рис. 11) і натиснемо OK (Transform selected

series). З’явиться графік ряду з виділеним трендом (див. рис. 12).

Рис. 12

Повернімось у вікно Transformation of variables. Перейдемо у закладку

Autocorrs і натиснемо кнопку Partial autocorrelations

(див. рис.

13).

Якщо значення змінної в різні моменту часу не пов’язані між собою, то

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

за межі вертикальних червоних ліній. Бачимо, що на нашому графіку (див.

рис. 14), за межі червоних ліній виходить тільки перший стовпчик, тобто

значення x_t залежить від xt-1, але майже не залежить від xt-2 , x_t-3, … .

Побудуємо модель для даного часового ряду. У закладці Autocorrs

натиснемо кнопку Autocorrelations

(див. рис.

13). Бачимо, що графік

кореляційної функції має типову спадну форму

(див. рис.

15). За

формою графіків часткової кореляційної функції та кореляційної функції

ми можемо зробити висновок, що поведінку ряду з виділеним трендом має

описувати модель авторегресії AR(1).

Повернімось у вікно Single Series ARIMA і в Arima model parameters у

полі p-Autoregressive вкажемо

(див. рис.

16). Натиснемо OK (Begin

parameter estimation).

Процес оцінки параметрів збіжний і з’являється вікно з результатами

аналізу

(див. рис.

17). Натиснувши у закладці Quick кнопку Summary:

Parameter estimates, отримаємо оцінки параметрів моделі. У нашому

випадку (див. рис. 18), при

1 tx − отримали коефіцієнт β ≈0.95. Це

означає, що наша модель має вигляд

=0.95

+ε

t−1

Рис. 17

Рис. 18

Щоб переконатися, що побудована модель правильна, треба

проаналізувати залишки. У вікні Single Series ARIMA Results перейдемо у

закладку Autocorrelation і у розділі Autocorrelation of residuals натиснемо

спочатку Autocorrelations, а потім Partial autocorrelations (див. рис. 19).

Бачимо, що на обох графіках

(див. рис.

20,

21) горизонтальні

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

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

Рис. 21

У закладці Distribution of Residuals

(див. рис. 22) можна подивитися

на узгодженість залишків з нормальним законом розподілу натиснувши

кнопки Histogram та Normal probability plot

(див. рис.

23,

24).

Бачимо, що нормальний розподіл досить добре описує залишки моделі.

Для того, щоб передбачити, як часовий ряд буде вести себе у

майбутньому, перейдемо у закладку Advanced і натиснемо кнопку Plot

series & forecasts

(див. рис. 25). При цьому на графіку (див. рис. 26)

червона лінія відповідатиме прогнозованим значенням, а зеленими лініями

буде виділено 90-відсотковий проміжок надійності передбачення.

Можна змінювати рівень надійності у полі Confidence level, вказувати

кількість передбачуваних значень у полі Number of cases, а також значення

з якого потрібно починати передбачення у віконці Start at case (див. рис.

25). Щоб отримати числові значення прогнозованих показників

натиснемо кнопку Forecast cases.

Часові ряди (перетворення

та інтервенції)

Відкриємо файл Series_g.sta, в якому містяться щомісячні дані про

кількість перевезень пасажирів авіакомпанією. Виконаємо Statistics Æ

Advanced Linear/ Nonlinear Models Æ Time Series/Forecasting. Побудуємо

графік ряду

(див. рис. 1). З графіка зрозуміло, що у часового ряду є

лінійний тренд та коливання, амплітуда котрих зростає із часом.

Спочатку позбавимося зростання амплітуди коливань. Для цього

прологарифмуємо ряд. Перейдемо у закладку Advanced натиснемо кнопку

Other transformations & plots. У панелі, що з’явиться, на закладці x=f(x),

відмітимо Natural log (x=ln(x)) і натиснемо OK (Transform selected series).

Отримаємо графік перетвореного часового ряду з постійною амплітудою

періодичної компоненти (див. рис. 2).

Далі позбудемося лінійного тренду. Для цього застосуємо до ряду

різницевий оператор із кроком 1. У закладці Difference, integrate відмітимо

Differencing (x=x-x(lag)) і вкажемо крок для різницевого оператора lag=1

(див. рис.

3). Натиснемо OK

(Transform selected series). Отримаємо

графік ряду з виділеним трендом (див. рис. 4).

Рис. 1

Тепер спробуємо задати ARIMA модель для даних, що отримали після

наведених вище перетворень. Повернемося у вікно Single series ARIMA та

вкажемо на закладці Quick в ARIMA model parameters p=0, P=0, q=1, Q=1.

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

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

більшими ніж три. Відмітимо у верхньому вікні наш початковий ряд і

вкажемо у Transform variable (series) prior to analysis перетворення, які ми

здійснили із рядом до аналізу: Natural log та Difference (див. рис. 6). N of

passes вказує скільки разів ми застосовували різницевий оператор.

Натискаємо OK (Begin Parameter estimation). Процес оцінки параметрів

збіжний. Далі, як і раніше, можемо дослідити кореляцію залишків та

пересвідчитись, що модель обрана вдало. Прогноз поведінки часового ряду

у майбутньому зображено на рисунку

Розглянемо приклад дослідження часового ряду з інтервенцією.

Відкриємо файл Director.sta, в якому містяться щомісячні дані про

кількість телефонних дзвінків, що надходили до компанії. Спочатку

дзвінки були безкоштовними, однак на 147 місяці компанія вирішила

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

Виконаємо Statistics Æ Advanced Linear/ Nonlinear Models Æ Time Series

/Forecasting Æ Interrupted time series analysis. Так як і раніше побудуємо

графік ряду (див. рис. 8).

Рис. 6

Оскільки на графіку видно різку зміну у значеннях спостережень -

інтервенцію, то перед побудовою моделі подивимось які типи інтервенції

можливі. Перейдемо на закладку Review Impact Patterns (див. рис. 9).

Змінюючи тип інтервенції у віконечку Type of intervention і натискаючи

кнопку Review types of impact patterns, отримуємо графічне зображення

різних типів інтервенції (див., наприклад, рис. 10). Для нашого випадку

підходить саме тип Abrupt, permanent, бо під час інтервенції

спостерігається різкий стрибок, а потім повільне зростання продовжуєтьсь

як і до моменту інтервенції.

Перейдемо на закладку Quick

(див. рис. 11) і у Specify times and types

of interventions вкажемо At case number - 147 і Type of intervention - Abrupt,

permanent. Вкажемо як параметри моделі p=0, P=0, q=0, Q=1. Для

виділення сезонної компоненти та тренду в Transform variable (series) prior

to analysis задамо Difference з кроком Lag = 1 і Lag=12 та в обох випадках

N of passes =1. Натискаємо OK (Begin parameter estimation). Процес оцінки

параметрів збіжний.

Далі аналіз результатів, залишків та адекватності моделі аналогічний до

попередніх прикладів.

Рис. 9

ття 12 Кластерний

аналіз

Відкриємо файл Cars.sta, який знаходиться у Examples/Datasets. У

даному файлі

(див. рис.

1) міститься інформація про характеристики

автомобілів різних виробників

(ціна, прискорення, надійність, легкість

управління, споживання пального). Завдання полягає у розбитті

автомобілів на групи за спільними ознаками.

Рис. 1

Для того, щоб мати можливість проаналізувати інформацію, треба

ввести шкалу порівняння для показників. Оскільки порівнювати ціну

автомобіля в доларах та прискорення в м/с² беззмістовно, то спочатку

змінні потрібно стандартизувати (тобто відняти середнє та поділити на

Виберемо Joining

(tree clustering). Натискаємо OK. У вікні, яке

з’явилось, на закладці Advanced (див. рис. 4), як Variables оберемо всі

змінні. Оскільки ми розбиваємо на групи виробників, тобто Cases, то у полі

Cluster виберемо Cases (rows). В полі Amalgamation (linkage rule) потрібно

обрати спосіб визначення відстані між кластерами.

Розбиття на кластери можна виконувати по-різному. Якщо ми знаємо,

що наші кластери мають форму довгих ланцюгів, то відстані між

кластерами зручно визначати, як відстані між найближчими точками

кластерів, тобто Single linkage. Якщо ж кластери мають форму плям, що

сконцентровані в кількох областях простору, то тоді доцільно визначити

відстань між кластерами, як відстань між найвіддаленішими точками

кластерів (Complete linkage), або їх центрами.

У нашому випадку виберемо Complete linkage. Як метод обчислення

відстані між двома точками в полі Distance measure виберемо Euclidean

distances - звичайну евклідову відстань. Натиснемо OK. З’явиться вікно, у

якому можна подивитися на результати аналізу (див. рис. 5).

Для того, щоб бачити, як відбувається розбиття на кластери, натиснемо

кнопку Horizontal hierarchical tree plot. Отримуємо дендрограму

ієрархічне дерево-графік, яке демонструє процес поступового об’єднання

досліджуваних об’єктів у кластери, відповідно до значень їхніх

характеристик (див. рис. 6).

Рис. 4

Для того, щоб подивитися, яким чином відбувався процес об’єднання в

кластери, можна у вікні Joining Results натиснути кнопку Amalgamation

Schedule і отримати табличку, в якій показано, як об’єкти спостереження

поступово об’єднувалися в кластери (див. рис. 7).

Рис. 7

Щоб дендрограма була вертикальною у вікні Joining Results слід

натиснути кнопку Vertical icicle plot.

Якщо хочемо бачити на дендрограмі відносні відстані, за яких точки

потраплять до кластерів, то у вікні Joining Results слід відмітити Scale tree

to dlink/dmax*100. Отримаємо вертикальну дендрограму з школою

відстаней об’єднання від 0 до 100 (див. рис. 8).

Вивчивши дендрограму, можемо зробити висновок, що кількість

кластерів складає 3 або 4.

Щоб мати змогу певним чином проінтерпретувати особливості груп, на

які були поділені виробники автомобілів, повернемось до вікна Сluster

Analysis і виконаємо K-means clustering

(див. рис.

3). У вікні, що

з’явиться

(див. рис. 9), перейдемо на закладку Advanced, виберемо всі

змінні і вкажемо Number of clusters - 3, а у віконці Cluster вкажемо Cases

(rows). Натиснемо OK.

У вікні k-Means Clustering Results, що з’явиться (див. рис. 10),

перейдемо на закладку Advanced і натиснемо кнопку Graph of means.

Рис. 10

Бачимо

(див. рис. 11), що до другого кластеру потрапили виробники

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

у керуванні та споживають мало пального. І навпаки, до третього кластеру

потрапили виробники недорогих ненадійних автомобілів, котрі

споживають багато пального, якими важко керувати, які, однак, мають

великий показник прискорення. До першого кластеру потрапили

виробники з середніми показниками.

Щоб побачити, які виробники потрапили до якого кластеру, натиснемо

кнопку Members of each cluster & distances. Результат зображено на

рисунку

142

Якщо ми хочемо кластеризувати не тільки виробників, а пари виробник-

характеристика (тобто кожна клітинка таблиці є точкою для кластеризації),

то у вікні Сluster Analysis слід вибрати Two way joining (див. рис. 3).

Обираємо для аналізу всі змінні і натискаємо OK. У вікні результатів

аналізу натискаємо кнопку Summary (див. рис. 13).

Рис. 13

На отриманому графіку (див. рис. 14) графічно зображено розбиття

на групи - кожен колір відповідає певному кластеру.

Two way joining методи застосовують у тих випадках, коли потрібно

розбити на групи не лише носіїв певних характеристик, але одночасно й

значення характеристик, які призводять до потрапляння у певний кластер.

Two way joining застосовують набагато рідше ніж інші методи кластерного

аналізу.

Заняття 13

Дискримінантний аналіз

Відкриємо файл Irisdat.sta, який знаходиться у Examples/Datasets. У

цьому файлі

(див. рис.

1) наведені дані Фішера про результати

вимірювань довжини і ширини чашолистків і пелюсток квітів іриса. Всього

було здійснено 150 вимірювань квіток іриса - по 50 для кожного з трьох

типів. Завдання полягає у тому, щоб маючи такі дані віднести новий ірис

до одного з цих трьох типів: SETOSA, VERSICOL, VIRGINIC. На відміну від

кластерного аналізу групи вже відомі, а дані, які використовує

дискримінантний аналіз, такі ж як і у кластерному аналізі - характеристики

представників цих груп. Завдання дискримінантного аналізу полягає у

встановлені за характеристиками нового об’єкта до якої групи він

належить.

Виконаємо Statistics Æ Multivariate Exploratory Techniques Æ

Discriminant Analysis

(див. рис. 2). У Variables виберемо як групуючу

змінну

(Grouping variable)

- змінну IRISTYPE, а як незалежні змінні

(Independent variables) - SEPALLEN, SEPALWD, PETALLEN, PETALWD. У

Codes for grouping variable обираємо All. Відмітимо також Advanced

options (stepwise analysis) (див. рис. 3).

Рис. 1