Введение

В данной лабораторной работе рассматриваются методы классификации и анализа данных в языке R.

В работе используются:

Для выполнения заданий используется встроенный набор данных iris. Этот набор содержит 150 наблюдений о цветках ириса. В нем есть четыре числовых признака и один столбец с видом цветка.

Перед запуском отчета необходимые пакеты нужно один раз установить в консоли RStudio:

install.packages("class")
install.packages("gmodels")
install.packages("caret")
install.packages("e1071")
install.packages("vegan")

Задание 2. Классификация k-ближайших соседей

Условие

Выполнить классификацию k-ближайших соседей с использованием функции knn() из пакета class на наборе данных iris.

Необходимо:

  • провести нормализацию данных;
  • разделить выборку на обучающую и тестовую;
  • оценить модель с использованием функции CrossTable() из пакета gmodels;
  • построить матрицу ошибок;
  • выполнить диагональную оценку качества прогноза.

Выполнение задания

# Подключаем пакет class.
# В нем находится функция knn(), которая выполняет классификацию
# методом k-ближайших соседей.

library(class)

# Подключаем пакет gmodels.
# Он нужен для функции CrossTable(), которая строит перекрестную таблицу
# реальных и предсказанных классов.

library(gmodels)

# Подключаем пакет caret.
# Он нужен для функции confusionMatrix(),
# которая строит матрицу ошибок и считает метрики качества.

library(caret)
# Загружаем встроенный набор данных iris.

data(iris)

# Смотрим структуру набора данных.
# str() показывает типы столбцов и общее устройство таблицы.

str(iris)
## 'data.frame':    150 obs. of  5 variables:
##  $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
##  $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
##  $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
##  $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
##  $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...
# summary() выводит краткую статистику.
# Для числовых столбцов показываются минимум, максимум, среднее, медиана и квартили.
# Для Species показывается количество объектов каждого вида.

summary(iris)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    :50  
##  versicolor:50  
##  virginica :50  
##                 
##                 
## 

В наборе iris содержится 150 наблюдений. Первые четыре столбца являются числовыми признаками, а столбец Species содержит класс цветка.

Нормализация данных

Метод KNN работает на основе расстояний между объектами. Поэтому признаки нужно привести к одному масштабу.

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

# Создаем функцию min-max нормализации.
# Она переводит значения признака в диапазон от 0 до 1.
#
# min(x) — минимальное значение признака.
# max(x) — максимальное значение признака.
# x - min(x) — сдвиг значений так, чтобы минимум стал равен 0.
# Деление на max(x) - min(x) приводит значения к масштабу от 0 до 1.

normalize <- function(x) {
  (x - min(x)) / (max(x) - min(x))
}
# Берем только числовые признаки.
# Столбец Species не нормализуем, потому что это класс, а не числовой признак.

iris_numeric <- iris[, 1:4]

# lapply() применяет функцию normalize к каждому столбцу.
# as.data.frame() превращает результат обратно в таблицу.

iris_norm <- as.data.frame(
  lapply(iris_numeric, normalize)
)

# Выводим первые строки нормализованных данных.

head(iris_norm)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1   0.22222222   0.6250000   0.06779661  0.04166667
## 2   0.16666667   0.4166667   0.06779661  0.04166667
## 3   0.11111111   0.5000000   0.05084746  0.04166667
## 4   0.08333333   0.4583333   0.08474576  0.04166667
## 5   0.19444444   0.6666667   0.06779661  0.04166667
## 6   0.30555556   0.7916667   0.11864407  0.12500000
# Проверим, как изменилась шкала признака Sepal.Length.
# До нормализации значения были в исходном диапазоне.

cat("До нормализации:\n")
## До нормализации:
cat("Минимум:", min(iris$Sepal.Length), "\n")
## Минимум: 4.3
cat("Максимум:", max(iris$Sepal.Length), "\n")
## Максимум: 7.9
# После нормализации минимум должен быть 0, максимум должен быть 1.

cat("После нормализации:\n")
## После нормализации:
cat("Минимум:", min(iris_norm$Sepal.Length), "\n")
## Минимум: 0
cat("Максимум:", max(iris_norm$Sepal.Length), "\n")
## Максимум: 1

Разделение на обучающую и тестовую выборки

# Фиксируем генератор случайных чисел.
# Это нужно, чтобы при повторном запуске отчета разбиение было одинаковым.

set.seed(123)

# Создаем случайные индексы строк для обучающей выборки.
# 70% данных пойдет на обучение модели.

train_index <- sample(
  1:nrow(iris_norm),
  size = 0.7 * nrow(iris_norm)
)

# Формируем обучающую выборку.
# На этих данных алгоритм будет искать ближайших соседей.

train_data <- iris_norm[train_index, ]

# Формируем тестовую выборку.
# На этих данных будем проверять качество классификации.

test_data <- iris_norm[-train_index, ]

# Сохраняем правильные классы для обучающей выборки.

train_labels <- iris$Species[train_index]

# Сохраняем правильные классы для тестовой выборки.

test_labels <- iris$Species[-train_index]

# Выводим размеры выборок.

cat("Обучающая выборка:", nrow(train_data), "строк\n")
## Обучающая выборка: 105 строк
cat("Тестовая выборка:", nrow(test_data), "строк\n")
## Тестовая выборка: 45 строк

Обучение и применение KNN

# Выбираем количество соседей k.
# В учебном примере можно взять корень из размера обучающей выборки.
# round() округляет результат до целого числа.

k_value <- round(sqrt(nrow(train_data)))

cat("Используем k =", k_value, "\n")
## Используем k = 10
# Выполняем классификацию методом KNN.
#
# train = train_data — обучающие данные.
# test = test_data — тестовые данные, для которых нужно предсказать класс.
# cl = train_labels — правильные классы обучающей выборки.
# k = k_value — количество ближайших соседей.

knn_result <- knn(
  train = train_data,
  test = test_data,
  cl = train_labels,
  k = k_value
)

# Выводим предсказанные классы.

knn_result
##  [1] setosa     setosa     setosa     setosa     setosa     setosa    
##  [7] setosa     setosa     setosa     setosa     setosa     setosa    
## [13] setosa     setosa     versicolor versicolor versicolor versicolor
## [19] versicolor versicolor versicolor versicolor versicolor versicolor
## [25] versicolor versicolor versicolor virginica  versicolor versicolor
## [31] versicolor versicolor virginica  virginica  virginica  virginica 
## [37] virginica  virginica  virginica  virginica  virginica  versicolor
## [43] virginica  virginica  virginica 
## Levels: setosa versicolor virginica

Оценка модели через CrossTable

# CrossTable() строит таблицу сравнения:
# x = реальные классы,
# y = предсказанные классы.
#
# Значения на главной диагонали — правильные ответы.
# Значения вне главной диагонали — ошибки классификации.
#
# prop.chisq = FALSE отключает лишнюю статистику хи-квадрат.

CrossTable(
  x = test_labels,
  y = knn_result,
  prop.chisq = FALSE
)
## 
##  
##    Cell Contents
## |-------------------------|
## |                       N |
## |           N / Row Total |
## |           N / Col Total |
## |         N / Table Total |
## |-------------------------|
## 
##  
## Total Observations in Table:  45 
## 
##  
##              | knn_result 
##  test_labels |     setosa | versicolor |  virginica |  Row Total | 
## -------------|------------|------------|------------|------------|
##       setosa |         14 |          0 |          0 |         14 | 
##              |      1.000 |      0.000 |      0.000 |      0.311 | 
##              |      1.000 |      0.000 |      0.000 |            | 
##              |      0.311 |      0.000 |      0.000 |            | 
## -------------|------------|------------|------------|------------|
##   versicolor |          0 |         17 |          1 |         18 | 
##              |      0.000 |      0.944 |      0.056 |      0.400 | 
##              |      0.000 |      0.944 |      0.077 |            | 
##              |      0.000 |      0.378 |      0.022 |            | 
## -------------|------------|------------|------------|------------|
##    virginica |          0 |          1 |         12 |         13 | 
##              |      0.000 |      0.077 |      0.923 |      0.289 | 
##              |      0.000 |      0.056 |      0.923 |            | 
##              |      0.000 |      0.022 |      0.267 |            | 
## -------------|------------|------------|------------|------------|
## Column Total |         14 |         18 |         13 |         45 | 
##              |      0.311 |      0.400 |      0.289 |            | 
## -------------|------------|------------|------------|------------|
## 
## 

Матрица ошибок и диагональная оценка качества

# confusionMatrix() строит матрицу ошибок.
# Первый аргумент — предсказанные классы.
# Второй аргумент — реальные классы.

knn_matrix <- confusionMatrix(
  knn_result,
  test_labels
)

knn_matrix
## Confusion Matrix and Statistics
## 
##             Reference
## Prediction   setosa versicolor virginica
##   setosa         14          0         0
##   versicolor      0         17         1
##   virginica       0          1        12
## 
## Overall Statistics
##                                           
##                Accuracy : 0.9556          
##                  95% CI : (0.8485, 0.9946)
##     No Information Rate : 0.4             
##     P-Value [Acc > NIR] : 2.842e-15       
##                                           
##                   Kappa : 0.9326          
##                                           
##  Mcnemar's Test P-Value : NA              
## 
## Statistics by Class:
## 
##                      Class: setosa Class: versicolor Class: virginica
## Sensitivity                 1.0000            0.9444           0.9231
## Specificity                 1.0000            0.9630           0.9688
## Pos Pred Value              1.0000            0.9444           0.9231
## Neg Pred Value              1.0000            0.9630           0.9688
## Prevalence                  0.3111            0.4000           0.2889
## Detection Rate              0.3111            0.3778           0.2667
## Detection Prevalence        0.3111            0.4000           0.2889
## Balanced Accuracy           1.0000            0.9537           0.9459
# Главная диагональ матрицы ошибок показывает количество правильных ответов
# по каждому классу.

diag_values <- diag(knn_matrix$table)

cat("Верно классифицировано по классам:\n")
## Верно классифицировано по классам:
diag_values
##     setosa versicolor  virginica 
##         14         17         12
# Диагональная оценка качества прогноза показывает,
# какая доля объектов каждого класса была классифицирована правильно.
#
# colSums(knn_matrix$table) считает общее количество объектов каждого реального класса.
# diag_values делится на это количество.

diag_quality <- diag_values / colSums(knn_matrix$table)

cat("Диагональная оценка качества по классам, %:\n")
## Диагональная оценка качества по классам, %:
round(diag_quality * 100, 1)
##     setosa versicolor  virginica 
##      100.0       94.4       92.3
# Общая точность модели показывает долю всех правильных ответов.

cat("Общая точность модели, %:\n")
## Общая точность модели, %:
round(knn_matrix$overall["Accuracy"] * 100, 2)
## Accuracy 
##    95.56

Подбор оптимального k

# Дополнительно проверим качество модели при разных значениях k.
# Берем только нечетные k от 1 до 19.
# Нечетные значения помогают избежать ничьей при голосовании соседей.

k_values <- seq(1, 20, by = 2)

# Для каждого k строим KNN и считаем точность.

accuracies <- sapply(k_values, function(k) {
  pred <- knn(
    train = train_data,
    test = test_data,
    cl = train_labels,
    k = k
  )

  mean(pred == test_labels)
})

# Собираем результаты в таблицу.

results_k <- data.frame(
  k = k_values,
  accuracy_percent = round(accuracies * 100, 1)
)

results_k
##     k accuracy_percent
## 1   1             95.6
## 2   3             95.6
## 3   5             95.6
## 4   7             95.6
## 5   9             97.8
## 6  11             97.8
## 7  13             97.8
## 8  15             97.8
## 9  17             97.8
## 10 19             97.8
# Строим график зависимости точности от k.
# Он помогает визуально выбрать подходящее значение k.

plot(
  k_values,
  accuracies * 100,
  type = "b",
  pch = 19,
  xlab = "k, число соседей",
  ylab = "Точность, %",
  main = "Выбор оптимального k для KNN",
  ylim = c(80, 100)
)

abline(
  h = max(accuracies * 100),
  lty = 2
)

Вывод по заданию 2

В задании была выполнена классификация методом k-ближайших соседей.

Перед применением KNN данные были нормализованы, потому что метод использует расстояния между объектами. Нормализация привела все признаки к диапазону от 0 до 1.

После этого данные были разделены на обучающую и тестовую выборки. С помощью функции knn() были получены предсказания для тестовых объектов.

Для оценки качества использовались CrossTable() и матрица ошибок. Значения на главной диагонали показали количество правильно классифицированных объектов. Диагональная оценка качества позволила определить точность отдельно по каждому классу.

В целом модель хорошо справилась с классификацией набора iris.

Задание 3. Метод опорных векторов

Условие

Рассмотреть пример реализации метода опорных векторов с использованием функции svm() из пакета e1071.

Необходимо:

  • построить линейный классификатор;
  • выполнить прогнозирование;
  • выполнить перекрестную проверку с делением исходной выборки на 10 равных частей (cross = 10).

Выполнение задания

# Подключаем пакет e1071.
# Он содержит функцию svm(), которая реализует метод опорных векторов.

library(e1071)

# caret нужен для матрицы ошибок.

library(caret)
# Фиксируем случайность для воспроизводимости результата.

set.seed(42)

# Делим данные на обучающую и тестовую выборки.
# 70% данных используется для обучения, 30% — для проверки.

svm_index <- sample(
  1:nrow(iris),
  size = 0.7 * nrow(iris)
)

train_svm <- iris[svm_index, ]
test_svm <- iris[-svm_index, ]

cat("Обучающая выборка:", nrow(train_svm), "строк\n")
## Обучающая выборка: 105 строк
cat("Тестовая выборка:", nrow(test_svm), "строк\n")
## Тестовая выборка: 45 строк

Построение линейного SVM-классификатора

# Строим SVM-модель.
#
# Species ~ . означает, что Species — целевая переменная,
# а все остальные столбцы используются как признаки.
#
# kernel = "linear" задает линейный классификатор.
# cost = 1 задает штраф за ошибки.
# scale = TRUE включает автоматическое масштабирование признаков.

svm_linear <- svm(
  Species ~ .,
  data = train_svm,
  kernel = "linear",
  cost = 1,
  scale = TRUE
)

svm_linear
## 
## Call:
## svm(formula = Species ~ ., data = train_svm, kernel = "linear", cost = 1, 
##     scale = TRUE)
## 
## 
## Parameters:
##    SVM-Type:  C-classification 
##  SVM-Kernel:  linear 
##        cost:  1 
## 
## Number of Support Vectors:  25
# nSV показывает количество опорных векторов по классам.
# sum() считает их общее количество.
#
# Опорные векторы — это точки, которые находятся ближе всего
# к границе между классами и сильнее всего влияют на модель.

cat("Число опорных векторов:", sum(svm_linear$nSV), "\n")
## Число опорных векторов: 25

Визуализация линейного классификатора

# Строим график SVM-модели.
# Для визуализации используются два признака:
# Petal.Length и Petal.Width.
#
# На графике видно, как модель делит пространство признаков
# на области классов.

plot(
  svm_linear,
  data = train_svm,
  formula = Petal.Width ~ Petal.Length,
  main = "SVM: линейный классификатор"
)

Предсказание и оценка модели

# predict() применяет обученную модель к тестовой выборке.
# В результате получаем предсказанные классы.

svm_prediction <- predict(
  svm_linear,
  test_svm
)

svm_prediction
##          7         11         12         19         23         25         28 
##     setosa     setosa     setosa     setosa     setosa     setosa     setosa 
##         31         37         39         45         46         51         52 
##     setosa     setosa     setosa     setosa     setosa versicolor versicolor 
##         56         59         62         63         66         67         70 
## versicolor versicolor versicolor versicolor versicolor versicolor versicolor 
##         75         78         79         82         91         95        101 
## versicolor  virginica versicolor versicolor versicolor versicolor  virginica 
##        102        106        108        112        113        116        117 
##  virginica  virginica  virginica  virginica  virginica  virginica  virginica 
##        119        125        127        133        134        137        141 
##  virginica  virginica  virginica  virginica versicolor  virginica  virginica 
##        143        147        148 
##  virginica  virginica  virginica 
## Levels: setosa versicolor virginica
# Строим матрицу ошибок для SVM.

svm_matrix <- confusionMatrix(
  svm_prediction,
  test_svm$Species
)

svm_matrix
## Confusion Matrix and Statistics
## 
##             Reference
## Prediction   setosa versicolor virginica
##   setosa         12          0         0
##   versicolor      0         14         1
##   virginica       0          1        17
## 
## Overall Statistics
##                                           
##                Accuracy : 0.9556          
##                  95% CI : (0.8485, 0.9946)
##     No Information Rate : 0.4             
##     P-Value [Acc > NIR] : 2.842e-15       
##                                           
##                   Kappa : 0.9324          
##                                           
##  Mcnemar's Test P-Value : NA              
## 
## Statistics by Class:
## 
##                      Class: setosa Class: versicolor Class: virginica
## Sensitivity                 1.0000            0.9333           0.9444
## Specificity                 1.0000            0.9667           0.9630
## Pos Pred Value              1.0000            0.9333           0.9444
## Neg Pred Value              1.0000            0.9667           0.9630
## Prevalence                  0.2667            0.3333           0.4000
## Detection Rate              0.2667            0.3111           0.3778
## Detection Prevalence        0.2667            0.3333           0.4000
## Balanced Accuracy           1.0000            0.9500           0.9537
# Выводим главную диагональ матрицы ошибок.
# Она показывает количество правильных ответов по каждому классу.

cat("Диагональ матрицы ошибок SVM:\n")
## Диагональ матрицы ошибок SVM:
diag(svm_matrix$table)
##     setosa versicolor  virginica 
##         12         14         17
# Выводим общую точность SVM-модели.

cat("Общая точность SVM, %:\n")
## Общая точность SVM, %:
round(svm_matrix$overall["Accuracy"] * 100, 2)
## Accuracy 
##    95.56

Перекрестная проверка cross = 10

# Выполняем кросс-валидацию.
#
# cross = 10 означает, что данные делятся на 10 частей.
# Модель 10 раз обучается и проверяется на разных частях данных.
#
# Это позволяет получить более надежную оценку качества модели.

svm_cross <- svm(
  Species ~ .,
  data = train_svm,
  kernel = "linear",
  cost = 1,
  scale = TRUE,
  cross = 10
)

cat("Точность по каждому из 10 разбиений:\n")
## Точность по каждому из 10 разбиений:
round(svm_cross$accuracies, 2)
##  [1] 100.00  90.91  90.00 100.00 100.00  90.91 100.00 100.00  90.00  90.91
# Считаем среднюю точность кросс-валидации.

cat("Средняя точность кросс-валидации:\n")
## Средняя точность кросс-валидации:
round(mean(svm_cross$accuracies), 2)
## [1] 95.27
# tot.acc — итоговая точность, которую возвращает svm() при cross > 0.

cat("Итоговая точность cross-validation:\n")
## Итоговая точность cross-validation:
round(svm_cross$tot.acc, 2)
## [1] 95.24

Подбор параметров SVM

# Дополнительно выполним подбор параметров.
# tune() перебирает разные значения cost и gamma
# и выбирает лучшую комбинацию по результатам кросс-валидации.

tune_result <- tune(
  svm,
  Species ~ .,
  data = train_svm,
  kernel = "radial",
  ranges = list(
    cost = c(0.1, 1, 10, 100),
    gamma = c(0.01, 0.1, 1)
  ),
  tunecontrol = tune.control(cross = 10)
)

cat("Лучшие параметры:\n")
## Лучшие параметры:
tune_result$best.parameters
##   cost gamma
## 3   10  0.01
cat("Ошибка лучшей модели:\n")
## Ошибка лучшей модели:
round(tune_result$best.performance * 100, 2)
## [1] 2.73
# Строим график ошибки при разных значениях параметров.
# Чем меньше ошибка, тем лучше модель.

plot(
  tune_result,
  main = "Ошибка SVM при разных cost и gamma"
)

# Используем лучшую модель, найденную tune().
# Проверяем ее на тестовой выборке.

best_svm <- tune_result$best.model

best_prediction <- predict(
  best_svm,
  test_svm
)

cat("Точность лучшей SVM-модели, %:\n")
## Точность лучшей SVM-модели, %:
round(mean(best_prediction == test_svm$Species) * 100, 2)
## [1] 97.78

Вывод по заданию 3

В задании был рассмотрен метод опорных векторов.

С помощью функции svm() была построена линейная модель классификации. Для обучения использовалась обучающая выборка, а качество проверялось на тестовой выборке.

Матрица ошибок показала, что SVM хорошо классифицирует объекты набора iris. Значения на диагонали матрицы ошибок соответствуют правильным предсказаниям.

Также была выполнена перекрестная проверка с cross = 10. Она позволяет оценить качество модели надежнее, чем одно случайное разбиение.

Дополнительно был выполнен подбор параметров через tune(). Это позволяет найти более подходящие значения параметров модели.

Задание 4. Метод главных компонент

Условие

Выполнить расчет главных компонент с использованием пакета vegan и функции rda().

Необходимо:

  • построить ординационную диаграмму методом PCA;
  • сделать выводы.

Выполнение задания

# Подключаем пакет vegan.
# В нем находится функция rda(), которую можно использовать для PCA.

library(vegan)

Подготовка данных

# Для PCA берем только числовые признаки.
# Species не используется, потому что это категориальная переменная.

iris_num <- iris[, 1:4]

head(iris_num)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1          5.1         3.5          1.4         0.2
## 2          4.9         3.0          1.4         0.2
## 3          4.7         3.2          1.3         0.2
## 4          4.6         3.1          1.5         0.2
## 5          5.0         3.6          1.4         0.2
## 6          5.4         3.9          1.7         0.4
# Выводим описательную статистику числовых признаков.

summary(iris_num)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500
# Считаем корреляционную матрицу.
# Она показывает связь между признаками.
#
# Значения, близкие к 1, говорят о сильной положительной связи.
# Значения, близкие к -1, говорят о сильной отрицательной связи.

round(cor(iris_num), 2)
##              Sepal.Length Sepal.Width Petal.Length Petal.Width
## Sepal.Length         1.00       -0.12         0.87        0.82
## Sepal.Width         -0.12        1.00        -0.43       -0.37
## Petal.Length         0.87       -0.43         1.00        0.96
## Petal.Width          0.82       -0.37         0.96        1.00

Расчет PCA через rda()

# rda() из пакета vegan может выполнять PCA.
#
# scale = TRUE обязательно используется для масштабирования признаков.
# Если признаки не масштабировать, признаки с большим разбросом
# могут сильнее влиять на результат.

pca_result <- rda(
  iris_num,
  scale = TRUE
)

pca_result
## 
## Call: rda(X = iris_num, scale = TRUE)
## 
##               Inertia Rank
## Total               4     
## Unconstrained       4    4
## 
## Inertia is correlations
## 
## Eigenvalues for unconstrained axes:
##    PC1    PC2    PC3    PC4 
## 2.9185 0.9140 0.1468 0.0207
# summary() показывает подробную информацию:
# собственные значения,
# долю объясняемой дисперсии,
# накопленную долю дисперсии.

summary(pca_result)
## 
## Call:
## rda(X = iris_num, scale = TRUE) 
## 
## Partitioning of correlations:
##               Inertia Proportion
## Total               4          1
## Unconstrained       4          1
## 
## Eigenvalues, and their contribution to the correlations 
## 
## Importance of components:
##                          PC1    PC2     PC3      PC4
## Eigenvalue            2.9185 0.9140 0.14676 0.020715
## Proportion Explained  0.7296 0.2285 0.03669 0.005179
## Cumulative Proportion 0.7296 0.9581 0.99482 1.000000

Объясняемая дисперсия

# eigenvals() возвращает собственные значения главных компонент.
# as.numeric() нужен, чтобы получить обычный числовой вектор.

eigen_values <- as.numeric(eigenvals(pca_result))

# Считаем долю объясняемой дисперсии в процентах.

variance_percent <- eigen_values / sum(eigen_values) * 100

# Считаем накопленную долю дисперсии.

cumulative_percent <- cumsum(variance_percent)

# Собираем результаты в таблицу.

variance_table <- data.frame(
  PC = paste0("PC", 1:length(eigen_values)),
  Eigenvalue = round(eigen_values, 3),
  Variance_percent = round(variance_percent, 2),
  Cumulative_percent = round(cumulative_percent, 2)
)

variance_table
##    PC Eigenvalue Variance_percent Cumulative_percent
## 1 PC1      2.918            72.96              72.96
## 2 PC2      0.914            22.85              95.81
## 3 PC3      0.147             3.67              99.48
## 4 PC4      0.021             0.52             100.00
# Строим график объясняемой дисперсии.
# Он показывает, сколько информации объясняет каждая компонента.

barplot(
  height = variance_percent,
  names.arg = paste0("PC", 1:length(variance_percent)),
  main = "Доля объясняемой дисперсии по компонентам",
  xlab = "Главная компонента",
  ylab = "Дисперсия, %",
  ylim = c(0, 100)
)

Ординационная диаграмма PCA

# Строим базовую ординационную диаграмму.
# На ней показано расположение объектов и направление признаков.

biplot(
  pca_result,
  display = c("sites", "species"),
  type = c("points", "text"),
  scaling = 2,
  main = "PCA biplot: объекты и переменные"
)

# Создаем цвета для разных видов ириса.

species_colors <- c(
  setosa = "red",
  versicolor = "green",
  virginica = "blue"
)

point_colors <- species_colors[as.character(iris$Species)]

# Получаем координаты объектов в пространстве главных компонент.

pca_scores <- scores(
  pca_result,
  display = "sites",
  scaling = 2
)

# Строим график первых двух главных компонент.

plot(
  pca_scores[, 1],
  pca_scores[, 2],
  col = point_colors,
  pch = 19,
  xlab = "PC1",
  ylab = "PC2",
  main = "PCA: ординационная диаграмма по видам"
)

legend(
  "topright",
  legend = names(species_colors),
  col = species_colors,
  pch = 19,
  title = "Вид"
)

Нагрузки признаков

# Нагрузки показывают вклад исходных признаков
# в формирование главных компонент.
#
# Чем больше значение по модулю,
# тем сильнее признак влияет на компоненту.

loadings_matrix <- scores(
  pca_result,
  display = "species",
  scaling = 0
)

round(loadings_matrix, 3)
##                 PC1    PC2
## Sepal.Length  0.521 -0.377
## Sepal.Width  -0.269 -0.923
## Petal.Length  0.580 -0.024
## Petal.Width   0.565 -0.067
## attr(,"const")
## [1] 4.940963
# Определим признаки, которые сильнее всего влияют на PC1 и PC2.

cat("Топ-признаки по PC1:\n")
## Топ-признаки по PC1:
sort(abs(loadings_matrix[, 1]), decreasing = TRUE)
## Petal.Length  Petal.Width Sepal.Length  Sepal.Width 
##    0.5804131    0.5648565    0.5210659    0.2693474
cat("Топ-признаки по PC2:\n")
## Топ-признаки по PC2:
sort(abs(loadings_matrix[, 2]), decreasing = TRUE)
##  Sepal.Width Sepal.Length  Petal.Width Petal.Length 
##   0.92329566   0.37741762   0.06694199   0.02449161

Вывод по заданию 4

В задании был выполнен расчет главных компонент с помощью функции rda() из пакета vegan.

Для анализа использовались четыре числовых признака набора iris. Столбец Species не использовался при расчете PCA, потому что он является категориальной переменной.

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

Ординационная диаграмма показывает, что вид setosa хорошо отделяется от остальных видов. Виды versicolor и virginica расположены ближе друг к другу, поэтому между ними возможно частичное пересечение.

По нагрузкам признаков видно, какие переменные сильнее всего влияют на главные компоненты. Признаки Petal.Length и Petal.Width сильно связаны с первой компонентой, а Sepal.Width заметно влияет на вторую компоненту.

Общий вывод

В лабораторной работе были рассмотрены методы классификации и анализа данных в языке R.

В задании 2 был применен метод k-ближайших соседей. Данные были нормализованы, разделены на обучающую и тестовую выборки, после чего была выполнена классификация. Качество модели было оценено через CrossTable(), матрицу ошибок и диагональную оценку качества прогноза.

В задании 3 был рассмотрен метод опорных векторов. Была построена линейная SVM-модель, выполнено прогнозирование и проведена перекрестная проверка с делением выборки на 10 частей.

В задании 4 был выполнен расчет главных компонент методом PCA с использованием функции rda() из пакета vegan. Были построены графики и сделаны выводы о структуре данных.

В результате работы были получены практические навыки классификации, оценки качества моделей, кросс-валидации и визуального анализа многомерных данных.