Кореляційний аналіз широко застосовується в дослідженнях геологічних даних. У цьому циклі статей ми розглянемо найбільш відомі методики розрахунку кореляційних залежностей: коефіцієнт кореляції r Пірсона, ранговий коефіцієнт кореляції p Спірмена і ранговий коефіцієнт кореляції τ Кендалла. Але почнемо з візуального контролю.
Перед візуалізацією та обчисленням коефіцієнтів кореляції давайте нагадаємо собі про те, що кореляція є мірою взаємозалежності показників. Мірою залежності одного показника від іншого є регресія, яку буде розглянуто в інших статтях. І ще - кореляція є мірою лінійної взаємозалежності між показниками. Відсутність кореляції ще не свідчить про відсутність інших типів взаємозалежності, наприклад - складної функціональної. Можна навести широко відомій графік із Вікіпедії, який показує різні типи взаємозалежності і відповідні значення коефіцієнтів кореляції:
 |
| Приклади взаємозалежності та кореляції |
Отже, перед обчисленням коефіцієнтів кореляції слід візуально визначити тип взаємозв’язків даних - побудувати діаграму розсіювання. Це можна зробити за допомогою функції plot:
>plot(granity$SiO2, granity$FeO, pch=20, cex=1.2, col=as.character(granity$color), main="Діаграма розсіювання", xlab="SiO2", ylab="FeO")#будуємо діаграму розсіювання
>plot(granity$SiO2, granity$FeO, pch=20, cex=1.2, col=as.character(granity$color), main="Діаграма розсіювання", xlab="SiO2", ylab="FeO")#будуємо діаграму розсіювання
 |
| Діаграма розсіювання |
Але у нашій таблиці 9 стовпців із числовими даними, які треба дослідити на наявність кореляційного зв’язку! Тобто нам треба побудувати (9*9)-9=72 графіків! "Чи не забагато доведеться клацати на клавіші?" - запитають у мене читачі, і будуть цілком праві. Для вирішення цієї проблеми в R існує безліч функцій, які будують комбіновані графіки із діаграмами розсіювання для всієї сукупності даних. Ми розглянемо функцію scatterplotMatrix базового пакету car:
>library(car)#завантажуємо пакет car
Loading required package: MASS
Loading required package: nnet#пакет car завантажує потрібні пакети
>scatterplotMatrix(granity[, 4:12], pch=20, cex=0.8)#будуємо графік із діаграмами розсіювання для стовпців 4-12 таблиці granity
>library(car)#завантажуємо пакет car
Loading required package: MASS
Loading required package: nnet#пакет car завантажує потрібні пакети
>scatterplotMatrix(granity[, 4:12], pch=20, cex=0.8)#будуємо графік із діаграмами розсіювання для стовпців 4-12 таблиці granity
(Я спеціально додав у команду побудови графіку параметри pch=20 та cex=0.8, щоб збільшити контрастність рисунку)
 |
| Діаграми розсіювання у функції scatterplotMatrix |
Вуаля! Функція доволі розумна - окрім діаграм розсіювання вона будує лінію регресії із довірчим інтервалом (червоні), графіки закону розподілу для кожного із компонентів (щільність вірогідності - сірі лінії у центральних віконцях) та лінію лінійної регресії (зелені). І це тільки базові параметри! Для більш детального ознайомлення я Вам раджу подивитись параметри цієї функції за допомогою команди ?scatterplotMatrix.
А тепер, давайте розглянемо цілком реальний випадок - коли даних настільки багато, що точки на діаграмі розсіювання наповзають одна на одну, зливаючись при цьому у суцільну пляму. Звісно, що при такому відображенні - майже нічого розгледіти не вдасться. Як можна викрутитись із подібної ситуації? Варіантів декілька. Наіпростіший - побудувати діаграму у графічний файл (наприклад - PNG) дуже великого розміру (маються на увазі сантиметри-міліметри, а не об’єм інформації). При цьому можна буде сподіватись, що точки на великому аркуші "розліпляться". Але це не шлях для професіоналів (А ми - справжні професіонали! Чи не так?). Все ж таки R розроблений науковцями-практиками для практичного використання. І з цією проблемою вони добре знайомі.
Існує багато методів відобразити діаграму розсіювання, яка характеризується величезною кількістю і сукупченістю точок. Один із найнаочніших методів - обчислити щільність точок на одиницю площі графіку і відобразити області з різною щільністю різними кольорами. Саме цей метод ми спробуємо, використавши для цього команду smoothScatter базового пакету graphics. І виключно для цього прикладу ми створемо дві вибірки, використавши генератор псевдовипадкових послідовностей та побудуємо діаграму розсіювання для цих вибірок:
>set.seed(1234)#встановлюємо первинне
#значення генератора псевдовипадкових послідовностей
>v1<- rnorm(10000, 50, 25)#генеруємо 10 тис. значень із середнім арифметичним 50 та середньоквадратичним відхиленням 25
>v2<- rnorm(8000, 50, 6)#генеруємо 8 тис. значень із середнім арифметичним 50 та середньоквадратичним відхиленням 6
>v3<- rnorm(2000, 40, 15)#генеруємо 2 тис. значень із середнім арифметичним 40 та середньоквадратичним відхиленням 15
>v4<- c(v2, v3)#об’єднуємо вибірки v2 та v3
>plot(v1,v4,pch=20)#будуємо діаграму розсіювання
>set.seed(1234)#встановлюємо первинне
#значення генератора псевдовипадкових послідовностей
>v1<- rnorm(10000, 50, 25)#генеруємо 10 тис. значень із середнім арифметичним 50 та середньоквадратичним відхиленням 25
>v2<- rnorm(8000, 50, 6)#генеруємо 8 тис. значень із середнім арифметичним 50 та середньоквадратичним відхиленням 6
>v3<- rnorm(2000, 40, 15)#генеруємо 2 тис. значень із середнім арифметичним 40 та середньоквадратичним відхиленням 15
>v4<- c(v2, v3)#об’єднуємо вибірки v2 та v3
>plot(v1,v4,pch=20)#будуємо діаграму розсіювання
 |
| Діаграма розсіювання псевдовипадкових вибірок |
Команду set.seed=1234 я використав для того, щоб задати первинний стан генератора псевдовипадкових послідовностей. Тепер Ви, задавши такий самий стан - отримаєте ідентичні зі мною вибірки (Довго шукав літературне джерело, але так і не знайшов - де бачив таку пропозицію: встановлювати різні первинні стани генератора псевдовипадкових послідовностей зручно викладачам при видачі практичних завдань студентам. Стан генератора при цьому відіграє роль "варіанта" завдання - і перевірити легко, і підробити не вдасться).
Тепер побудуємо графік щільності точок на діаграмі розсіювання:
>smoothScatter(v1, v4, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, nrpoints=500, main="Щільність точок на діаграмі розсіювання\nдля псевдовипадкових вибірок")
>smoothScatter(v1, v4, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, nrpoints=500, main="Щільність точок на діаграмі розсіювання\nдля псевдовипадкових вибірок")
 |
| Графік щільності точок на діаграмі розсіювання для псевдовипадкових вибірок |
Розберемо параметри команди smoothScatter: v1 та v4 - вибірки, за якими потрібно будувати графік; colramp - присвоювання кольорів (тут ми використали функцію colorRampPalette для змішування кольорів); pch задає тип відображення об’єктів на графіку. Зверніть увагу на параметр nrpoints - він задає кількість точок, які треба відобразити на графіку, починаючи від області з найменшою щільністю. За умовчанням цей параметр дорівнює 100 точкам (nrpoints=100). Щоб вилучити обмеження на кількість точок (тобто винести на діаграму всі точки) - треба цьому параметру присвоїти значення нескінченності: nrpoints=Inf. Але будьте обережні: при високій щільності точок - ви знов отримаєте чорну пляму всередені діаграми, і зафарбовані ділянки там, де щільність точок дозволить прорисовати кольори. І ще я хочу звернути Вашу увагу на головний підпис графіка: між словами "розсіювання" та "для" стоїть команда "\n", яка примусово починає новий рядок у підписі.
А тепер спробуємо застосувати функцію smoothScatter до наших реальних даних:
>smoothScatter(granity$SiO2, granity$Fe2O3, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, xlab="SiO2", ylab="Fe2O3", main="Щільність точок на діаграмі розсіювання")
>smoothScatter(granity$SiO2, granity$Fe2O3, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, xlab="SiO2", ylab="Fe2O3", main="Щільність точок на діаграмі розсіювання")
 |
| Графік щільності точок на діаграмі розсіювання |
І, наостаннє, - не забувайте про проблему масштабування. В принципі - побудовані вище діаграми розсіювання не передають справжнє відношення між показниками, оскільки використовують різні масштаби для осей абсцис та ординат. Через це відбувається розтягування однієї осі та стиснення іншої. Якщо ви хочете показати на діаграмі справжнє співвідношення між компонентами, вона повинна будуватись із застосуванням обмежень по осях - мінімальне значення осей повинно бути трохи менше за найменше значення з обох вибірок, а максимальне - трохи більше за найбільше значення:
>summary(granity$SiO2)#перевіряємо показники компоненту SiO2
Min. 1st Qu. Median Mean 3rd Qu. Max.
31.70 54.34 66.49 63.00 72.27 77.23
>summary(granity$Fe2O3)#перевіряємо показники компоненту Fe2O3
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
0.050 0.940 1.500 2.435 3.180 14.570 3
>smoothScatter(granity$SiO2, granity$Fe2O3, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, xlab="SiO2", ylab="Fe2O3", main="Щільність точок на діаграмі розсіювання\nіз масштабуванням", xlim=c(0,80), ylim=c(0,80))
>summary(granity$SiO2)#перевіряємо показники компоненту SiO2
Min. 1st Qu. Median Mean 3rd Qu. Max.
31.70 54.34 66.49 63.00 72.27 77.23
>summary(granity$Fe2O3)#перевіряємо показники компоненту Fe2O3
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
0.050 0.940 1.500 2.435 3.180 14.570 3
>smoothScatter(granity$SiO2, granity$Fe2O3, colramp=colorRampPalette(c("white", "darkgreen", "green", "blue", "magenta", "red", "yellow")), pch=20, xlab="SiO2", ylab="Fe2O3", main="Щільність точок на діаграмі розсіювання\nіз масштабуванням", xlim=c(0,80), ylim=c(0,80))
 |
| Графік щільності точок на діаграмі розсіювання із масштабуванням осей |
Рисунок не такий привабливий, як попередні. Але на ньому обидві вісі побудовані в однаковому масштабі, і дивлячись на цей графік ми можемо візуально "прикинути" - якого компоненту більше, на скільки більше, і в яких пропорціях компоненти містяться у досліджених породах.
У наступних статтях ми розглянемо "числові" методи кореляційного аналізу даних у програмному середовищі R, та проблеми статистичної значущості обчислених коефіцієнтів кореляції.
Немає коментарів:
Дописати коментар