Факторный анализ и анализ главных компонент

Проклятие размерности (curse of dimensionality)

С данными могут случиться две проблемы: либо их слишком мало, либо их слищком много. С первой проблемой многие так или иначе сталкивались — во-первых, не достигается требуемая статистическая мощность анализа (мало наблюдений), во-вторых, упускаем какие-либо закономерности (мало переменных). А в чем может быть проблема, если данных много?

Геометрическая вероятность

В единичный квадрат кидают точку. Какова вероятность, что она попадет во вписанный круг?

\[ p = \frac{S_{\mathrm{circ}}}{S_{\mathrm{sq}}} \\ p = \frac{\pi R^2}{1} = \pi \left( \frac{1}{2} \right)^2 = \frac{\pi}{4} \approx 0.79 \]

Окей, это для \(\mathbb{R}^2\), то есть для пространства размерности два (оно же — плоскость). Теперь рассмотрим ситуацию для \(\mathbb{R}^3\), то есть для пространства размерности три.

\[ p = \frac{V_{\mathrm{ball}}}{V_{\mathrm{cube}}} \\ p = \frac{\frac{4}{3} \pi R^3}{1} = \frac{4}{3} \pi \left( \frac{1}{2} \right)^3 = \frac{4 \pi}{24} \approx 0.52 \]

Хм, как-то маловато…

В общем случае, для \(\mathbb{R}^k\) объем шара равен

\[ k = 2n, \; V = \frac{\pi^n}{n!}R^{2n} \\ k = 2n + 1, \; V = \frac{2\cdot(2\pi)^n}{(2n+1)!!}R^{2n+1} \]

Можно аналитически доказать, что при \(k \rightarrow 0\): \(V \rightarrow 0\).

---

То есть, говоря более осязаемо, места вне шара становится больше. Набирая переменных, мы делаем расстояния между точками больше, из-за объем выборки должен расти экспоненциально, чтобы сохранялась адекватная точность оценки параметров. Это и называется проклятием размерности.

Что делать? Снижать размерность.

Факторный анализ vs анализ главных компонент

Так сложилось, что под термином «факторный анализ» часто смешиваются две вещи: анализ главных компонент и собственно факторный анализ. Кроме того, собственно факторный анализ также подразделяется на эксплораторный — когда мы ищем факторы — и конфирматорный — когда у нас есть гипотеза о структуре факторов, и нас необходимо ее подтвердить. Последний вариант встречается реже, поэтому мы будем рассматривать конфирматорный вариант факторного анализа.

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

Задачи факторного анализа и анализа главных компонент

1. Сокращение числа переменных
2. Измерение неизмеримого (построение новых обобзеных показателей)
3. Наглядное представление многомерных наблюдений
4. Описание структуры взаимных связей между переменными
5. Преодоление мультиколлинеарности (в регрессионном анализе)
6. Заполнение пропущенных значений (при работе с разряженными матрицами)

и т.д.

Анализ главных компонент (principal component analysis, PCA)

Математическая модель анализа главных компонент

Математические детали

Рассмотрим случайный вектор (матрицу) \(\begin{pmatrix} \boldsymbol{X}_1 & \boldsymbol{X}_2 & \dots & \boldsymbol{X}_k \end{pmatrix}\) (\(\boldsymbol{X}_i\) — некоторый столбец [числовых] данных). Мы хотим найти такую линейную комбинацию наших данных, у которой дисперсия максимальна. При это мы хотим домножать матрицу данных на вектор коэффициентов единичной длины, чтобы искусственно не раздувать дисперсию.

Формально, \[ \boldsymbol{Y}_1 = a_{11} \boldsymbol{X}_1 + a_{12} \boldsymbol{X}_2 + \dots + a_{1k} \boldsymbol{X}_k \\ \mathrm{var}(\boldsymbol{Y}_1) \rightarrow \max, \; \boldsymbol{a}_1 \boldsymbol{a}_1^\mathrm{T}= 1, \; \boldsymbol{a}_1 = (a_{11}, a_{12}, \dots, a_{1k}) \]

Далее, находим вторую линейную комбинацию наших данных с максимальной дисперсией, независимую (некоррелированную) с первой. Требования к вектору коэффициентов оставляем те же.

\[ \boldsymbol{Y}_2 = a_{21} \boldsymbol{X}_1 + a_{22} \boldsymbol{X}_2 + \dots + a_{2k} \boldsymbol{X}_k \\ \mathrm{var}(\boldsymbol{Y}_2) \rightarrow \max, \; \boldsymbol{a}_2 \boldsymbol{a}_2^\mathrm{T}= 1, \; \boldsymbol{a}_2 = (a_{21}, a_{22}, \dots, a_{2k}) \\ \mathrm{corr}(\boldsymbol{Y}_2, \boldsymbol{Y}_1) = 0 \]

Далее, находим третью линейную комбинацию, аналогичную предыдущим и некоррелированную с ними. \[ \boldsymbol{Y}_3 = a_{31} \boldsymbol{X}_1 + a_{32} \boldsymbol{X}_2 + \dots + a_{3k} \boldsymbol{X}_k \\ \mathrm{var}(\boldsymbol{Y}_3) \rightarrow \max, \; \boldsymbol{a}_3 \boldsymbol{a}_3^\mathrm{T}= 1, \; \boldsymbol{a}_3 = (a_{31}, a_{32}, \dots, a_{3k}) \\ \mathrm{corr}(\boldsymbol{Y}_3, \boldsymbol{Y}_1) = 0, \; \mathrm{corr}(\boldsymbol{Y}_3, \boldsymbol{Y}_2) = 0 \]

Так как у нас всего \(k\) переменных, то мы можем найти \(k\) таких линейных комбинаций.

В общем виде, они будут выглядеть следующим образом:

\[ \boldsymbol{Y}_k = a_{k1} \boldsymbol{X}_1 + a_{k2} \boldsymbol{X}_2 + \dots + a_{kk} \boldsymbol{X}_k \\ \mathrm{var}(\boldsymbol{Y}_k) \rightarrow \max, \; \boldsymbol{a}_k \boldsymbol{a}_k^\mathrm{T}= 1, \; \boldsymbol{a}_3 = (a_{31}, a_{32}, \dots, a_{3k}) \\ \forall \, 0 < i \leq k \; \mathrm{corr}(Y_k, Y_i) = 0 \]

Полученные \(\boldsymbol{Y}_i\) и будут искомые главные компоненты — наши новые оси, с помощью которых мы будем смотреть на данные и описывать их.

Если мы в качестве информативности компоненты используем её дисперсию, то можно сказать, что мы ищем наиболее информативные линейные комбинации. То есть, мы ищем некоторую «правильную систему координат».

Чтобы лучше понять, что все это значит, посмотрим на картинку.

Рассмотрим самый простой случай. Пусть у нас есть две переменные \((\boldsymbol{X}_1, \boldsymbol{X}_2)\) и мы подбираем их линейную комбинацию с наибольшей дисперсией. Графически это то же самое, что провести прямую \(\boldsymbol{Y}_1\), то есть такую прямую, по которой наши данные будут максимально «растянуты».

Далее ищем вторую линейную комбинацию. Прицнип тот же: необходима прямая, по которой данные макисмально растянуты. Однако эта прямая должна быть перпендикулярная \(\boldsymbol{Y}_1\), чтобы выполнялось условие \(\mathrm{corr}(\boldsymbol{Y}_1, \boldsymbol{Y}_2) = 0\) (в двумерном пространстве отсутствие корреляции эквивалентно ортогональности).

Получаем что-то такое. Красные оси и есть наши главные компоненты. Их столько же, сколько изначально было переменных. Теперь мы можем рассматривать наши наблюдения относительно новой системы координат \(Y_1 O Y_2\).

Сокращение размерности признакового пространства

Нуок, новые координаты мы нашли, но чё-то переменных, собственно, меньше не стало. А мы вапще-то хотели, чтобы осей стало поменьбше.

Идея состоит в следующем. Представим лаваш в трехмерном пространстве. Координатные оси пройдут через оси эллипсоида. Если мы посмотрим на лаваш как на признаковое пространство, то обнаружим, что большая часть вариативности данных распределена по осям \(x\) и \(y\), а вариативность по оси \(z\) можно отбросить, поскольку она мала и, по сути, примерно ошибка.

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

Реализация PCA в R

Грузим нужные пакеты:

library(tidyverse)

## ── Attaching packages ────────────────────────────────────────────────── tidyverse 1.2.1 ──

## ✔ ggplot2 3.3.0     ✔ purrr   0.3.2
## ✔ tibble  2.1.3     ✔ dplyr   0.8.3
## ✔ tidyr   1.0.2     ✔ stringr 1.4.0
## ✔ readr   1.3.1     ✔ forcats 0.4.0

## ── Conflicts ───────────────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()

library(corrplot)

## corrplot 0.84 loaded

theme_set(theme_bw())

Грузим данные:

pizza <- read_csv('https://raw.githubusercontent.com/angelgardt/hseuxlab-andan/master/Pizza.csv')

## Parsed with column specification:
## cols(
##   brand = col_character(),
##   id = col_double(),
##   mois = col_double(),
##   prot = col_double(),
##   fat = col_double(),
##   ash = col_double(),
##   sodium = col_double(),
##   carb = col_double(),
##   cal = col_double()
## )

str(pizza)

## Classes 'spec_tbl_df', 'tbl_df', 'tbl' and 'data.frame': 300 obs. of  9 variables:
##  $ brand : chr  "A" "A" "A" "A" ...
##  $ id    : num  14069 14053 14025 14016 14005 ...
##  $ mois  : num  27.8 28.5 28.4 30.6 30.5 ...
##  $ prot  : num  21.4 21.3 20 20.1 21.3 ...
##  $ fat   : num  44.9 43.9 45.8 43.1 41.6 ...
##  $ ash   : num  5.11 5.34 5.08 4.79 4.82 4.92 4.71 5.28 5.02 5.16 ...
##  $ sodium: num  1.77 1.79 1.63 1.61 1.64 1.65 1.58 1.75 1.71 1.66 ...
##  $ carb  : num  0.77 1.02 0.8 1.38 1.76 1.4 1.77 2.95 1.18 0.64 ...
##  $ cal   : num  4.93 4.84 4.95 4.74 4.67 4.67 4.63 4.72 4.93 4.95 ...
##  - attr(*, "spec")=
##   .. cols(
##   ..   brand = col_character(),
##   ..   id = col_double(),
##   ..   mois = col_double(),
##   ..   prot = col_double(),
##   ..   fat = col_double(),
##   ..   ash = col_double(),
##   ..   sodium = col_double(),
##   ..   carb = col_double(),
##   ..   cal = col_double()
##   .. )

Данные содержат следующие переменные:

• brand — Pizza brand (class label)
• id — Sample analysed
• mois — Amount of water per 100 grams in the sample
• prot — Amount of protein per 100 grams in the sample
• fat — Amount of fat per 100 grams in the sample
• ash — Amount of ash per 100 grams in the sample
• sodium — Amount of sodium per 100 grams in the sample
• carb — Amount of carbohydrates per 100 grams in the sample
• cal — Amount of calories per 100 grams in the sample

Поправляем данные:

pizza %>% mutate(brand = as_factor(brand),
             id = as_factor(id)) %>% 
  select(-brand, -id) -> pizza2 # убираем номинативные переменные

Смотрим описательные статистики:

summary(pizza2)

##       mois            prot            fat             ash       
##  Min.   :25.00   Min.   : 6.98   Min.   : 4.38   Min.   :1.170  
##  1st Qu.:30.90   1st Qu.: 8.06   1st Qu.:14.77   1st Qu.:1.450  
##  Median :43.30   Median :10.44   Median :17.14   Median :2.225  
##  Mean   :40.90   Mean   :13.37   Mean   :20.23   Mean   :2.633  
##  3rd Qu.:49.12   3rd Qu.:20.02   3rd Qu.:21.43   3rd Qu.:3.592  
##  Max.   :57.22   Max.   :28.48   Max.   :47.20   Max.   :5.430  
##      sodium            carb             cal       
##  Min.   :0.2500   Min.   : 0.510   Min.   :2.180  
##  1st Qu.:0.4500   1st Qu.: 3.467   1st Qu.:2.910  
##  Median :0.4900   Median :23.245   Median :3.215  
##  Mean   :0.6694   Mean   :22.865   Mean   :3.271  
##  3rd Qu.:0.7025   3rd Qu.:41.337   3rd Qu.:3.520  
##  Max.   :1.7900   Max.   :48.640   Max.   :5.080

Отмечаем разницу единиц измерения переменных — необходима стандартизация.

Смотрим разброс по переменным:

apply(pizza2, 2, sd)

##       mois       prot        fat        ash     sodium       carb 
##  9.5529866  6.4343924  8.9756583  1.2697235  0.3703577 18.0297225 
##        cal 
##  0.6200343

Разброс, ожидаемо, тоже различается — точно необходимо стандартизация.

Смотрим корреляции:

round(cor(pizza2), 2)

##         mois  prot   fat   ash sodium  carb   cal
## mois    1.00  0.36 -0.17  0.27  -0.10 -0.59 -0.76
## prot    0.36  1.00  0.50  0.82   0.43 -0.85  0.07
## fat    -0.17  0.50  1.00  0.79   0.93 -0.64  0.76
## ash     0.27  0.82  0.79  1.00   0.81 -0.90  0.33
## sodium -0.10  0.43  0.93  0.81   1.00 -0.62  0.67
## carb   -0.59 -0.85 -0.64 -0.90  -0.62  1.00 -0.02
## cal    -0.76  0.07  0.76  0.33   0.67 -0.02  1.00

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

Смотрим корреляции:

corrplot(cor(pizza2), diag = FALSE)

И еще разок:

GGally::ggpairs(pizza2)

## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2

Приступим собственно к PCA. Для этого есть команда prcomp.

pca <- prcomp(pizza2, scale. = TRUE) # обязательно стандартизовать переменные
pca

## Standard deviations (1, .., p=7):
## [1] 2.042494038 1.513425713 0.643865158 0.308503205 0.166364113 0.018374149
## [7] 0.003085252
## 
## Rotation (n x k) = (7 x 7):
##                PC1        PC2         PC3        PC4          PC5
## mois    0.06470937 -0.6282759  0.42166894  0.2207216  0.006470293
## prot    0.37876090 -0.2697067 -0.74602744  0.0105932  0.387982788
## fat     0.44666592  0.2343791  0.19930871  0.5070422 -0.173367634
## ash     0.47188953 -0.1109904 -0.05627269 -0.5523985 -0.670885701
## sodium  0.43570289  0.2016617  0.45516887 -0.4462769  0.602614079
## carb   -0.42491371  0.3203121 -0.05223651 -0.3343395 -0.007436899
## cal     0.24448730  0.5674576 -0.11331559  0.2792632 -0.078003175
##                  PC6           PC7
## mois   -0.4464499018 -0.4185690354
## prot    0.0001715203 -0.2767646428
## fat     0.5254028685 -0.3776715255
## ash    -0.0588609281 -0.0560214003
## sodium -0.0031309852  0.0005243238
## carb    0.0005088535 -0.7760679112
## cal    -0.7219138527 -0.0120598098

Собственно в аутпуте есть все, что нам нужно — корелляционная матрица главных компонент и переменных и дисперсии главных компонент. Но можно дисперсии представить в более удобоваримом варианте.

summary(pca)

## Importance of components:
##                          PC1    PC2     PC3    PC4     PC5     PC6
## Standard deviation     2.042 1.5134 0.64387 0.3085 0.16636 0.01837
## Proportion of Variance 0.596 0.3272 0.05922 0.0136 0.00395 0.00005
## Cumulative Proportion  0.596 0.9232 0.98240 0.9960 0.99995 1.00000
##                             PC7
## Standard deviation     0.003085
## Proportion of Variance 0.000000
## Cumulative Proportion  1.000000

Итак, нам будут нужны стандартные отклонения и cumulative proportion [of variance]. Как теперь отбросить компоненты?

Логический вариант следующий. Так как мы сдантартизовали переменные, то их дисперсии (стандартные отклонения) равны единице. Если дисперсия (стандартное) главной компоненты меньше единицы, значит она содержит информации меньше, чем исходная переменная, и её можно отбросить. То есть, в данно случае мы отбрасываем все компоненты, кроме \(\mathrm{PC1}\) и \(\mathrm{PC2}\).

Более аккуратный вариант, основанный на численных экспериментах, диктует нас правило, что нужно оставить все главные компоненты, дисперсия которых больше 0.8 (то есть стандартное отклонение больше 0.89).

Можно сделать визуализации:

ggplot(NULL, aes(names(summary(pca)$importance[2,]), summary(pca)$importance[2,])) +
  geom_point(size = 2) +
  geom_line(aes(group = 1)) +
  xlab('Principal components') + ylab('Variance Proportion')

ggplot(NULL, aes(names(summary(pca)$importance[3,]), summary(pca)$importance[3,])) +
  geom_point(size = 2) +
  geom_line(aes(group = 1)) +
  xlab('Principal components') + ylab('Cumulative Proportion')

Окей, мы сократили размерность. Как это использовать дальше? Нам нужны значения главных компонент. Их можно вытащить следующим образом:

pca$x

##               PC1        PC2           PC3           PC4          PC5
##   [1,]  5.0019853  2.6747462  3.926121e-02 -0.1438651878  0.283541584
##   [2,]  5.0153755  2.5250755  9.689453e-02 -0.3533598630  0.214909803
##   [3,]  4.7974235  2.6692404  7.534917e-02  0.1076180749 -0.034979111
##   [4,]  4.4620879  2.2812177  1.200328e-01  0.0539379991  0.174209451
##   [5,]  4.4644332  2.1555508  7.362777e-04 -0.1169702238  0.312504371
##   [6,]  4.4972855  2.1643567  1.747249e-01 -0.1152762515  0.200465801
##   [7,]  4.3081595  2.0536120 -7.545342e-06 -0.0164032892  0.265809424
##   [8,]  4.7499153  2.3491608  1.040862e-02 -0.4917731815  0.249287698
##   [9,]  4.8465215  2.6767466  1.012320e-01 -0.0193430903  0.171844126
##  [10,]  4.9082055  2.6540786 -6.868420e-02 -0.0006299387  0.053026307
##  [11,]  4.8329208  2.7320436 -8.834757e-02  0.1238664559  0.143665050
##  [12,]  4.8832096  2.7379271 -1.020999e-01 -0.0496339202  0.076776342
##  [13,]  4.6030631  2.5749513  1.607627e-03  0.2016839982 -0.051186621
##  [14,]  4.4141572  2.2878758  2.265809e-01 -0.0694073743  0.213881910
##  [15,]  4.4593259  2.1293271  4.750100e-02 -0.0784327295  0.183766393
##  [16,]  4.7455057  2.4630211  5.793408e-01 -0.2531866344  0.006169928
##  [17,]  4.4377161  1.8253153  2.766745e-01 -0.5168887260  0.161756594
##  [18,]  4.1950087  1.7984783  5.894449e-01 -0.2350983578 -0.039784548
##  [19,]  4.4893709  1.8829041  1.005629e-01 -0.4529566363 -0.072912335
##  [20,]  4.5998760  2.2313444  3.630431e-01 -0.3681811124  0.043173595
##  [21,]  4.8491661  2.9613955  1.032518e-01  0.2386910947 -0.004553835
##  [22,]  4.8146080  2.7820639  1.284386e-01  0.0405004270  0.138662806
##  [23,]  4.5108046  2.1588240  9.360943e-02 -0.1108467627  0.167457688
##  [24,]  4.2468605  2.1208991  1.858281e-01  0.3115869348  0.031688823
##  [25,]  4.4822615  2.0146307  5.759251e-01 -0.3409221361  0.015285109
##  [26,]  4.5340226  2.3761157  1.716958e-01  0.1446156869  0.118561799
##  [27,]  4.6263878  2.6054504  9.551799e-02  0.2072413215  0.123765355
##  [28,]  4.4613476  2.4298758  2.980755e-01  0.0648774943  0.043400727
##  [29,]  4.7424245  2.3570308  5.515856e-02 -0.1342822474  0.214202048
##  [30,]  1.8580800 -0.4850685  1.070249e+00  0.1276829184  0.034395245
##  [31,]  1.3101550 -1.1243616  9.127230e-01  0.2482539430  0.091073518
##  [32,]  1.7627285 -0.3056980  9.651271e-01  0.3826054752 -0.181362730
##  [33,]  1.6686052 -0.6503616  9.638492e-01  0.2225623629 -0.070703834
##  [34,]  1.5906089 -0.7258585  1.040190e+00  0.1159253359 -0.125675245
##  [35,]  1.7492185 -0.5240380  1.077354e+00  0.2452598341 -0.016154908
##  [36,]  1.7946998 -0.5735792  1.148347e+00  0.3089860545 -0.128490798
##  [37,]  1.2839956 -0.5866064  9.455694e-01  0.8431973574 -0.191863017
##  [38,]  0.9837179 -1.2976022  9.507788e-01  0.4408490449 -0.083543278
##  [39,]  1.5930404 -1.0780442  1.009201e+00 -0.2623823485  0.002470653
##  [40,]  1.6641408 -0.6360677  9.210032e-01  0.4210841573 -0.072082228
##  [41,]  1.6817448 -0.5023713  9.525960e-01  0.4158632671 -0.067709575
##  [42,]  1.7950103 -0.7891348  1.100283e+00 -0.0459552093  0.075272188
##  [43,]  1.7476960 -0.5125858  1.083379e+00  0.1664065441 -0.037018972
##  [44,]  1.0155957 -1.5499123  9.730668e-01  0.0976593721  0.154956456
##  [45,]  0.9593045 -1.6951569  8.414654e-01  0.0139594922  0.102425411
##  [46,]  1.0586805 -1.3616208  8.835314e-01  0.3859847407 -0.016534442
##  [47,]  1.6910284 -0.5107555  8.761171e-01  0.3933494419 -0.164888117
##  [48,]  1.8282608 -0.6252785  1.173273e+00  0.0750532084 -0.058673320
##  [49,]  1.5272432 -1.0874582  1.128503e+00 -0.2340756145  0.252895338
##  [50,]  1.6199850 -0.5797333  9.805356e-01  0.3100090549 -0.096610537
##  [51,]  1.7868292 -0.4510867  1.086624e+00  0.2106960693 -0.055592226
##  [52,]  1.7247135 -0.4180485  1.066904e+00  0.2660729588 -0.055407855
##  [53,]  1.7893329 -0.7429928  1.098574e+00 -0.0184825066  0.040787921
##  [54,]  1.7728750 -0.7821557  1.074858e+00 -0.0055463931  0.003666647
##  [55,]  1.2183717 -1.3339514  9.916531e-01  0.1589656669  0.057059326
##  [56,]  1.5412395 -0.9729048  1.044162e+00  0.1557772762 -0.070412522
##  [57,]  1.7061501 -0.8692735  1.067252e+00  0.1326282275  0.009950617
##  [58,]  1.7548225 -0.7927958  1.100360e+00  0.1562259906 -0.013591341
##  [59,]  1.2090193 -0.8094737  9.102128e-01  0.8014833039 -0.112046133
##  [60,]  1.2797590 -0.8030935  9.615940e-01  0.6419378690 -0.053466223
##  [61,]  1.3710443 -1.7982189 -1.243080e+00  0.4021198737  0.088162837
##  [62,]  1.1510979 -2.2500358 -1.654857e+00  0.1086459815  0.178447193
##  [63,]  1.1426014 -1.8234364 -1.161053e+00  0.4770536958  0.120913314
##  [64,]  1.2216241 -1.9344203 -1.138455e+00  0.1228446757  0.179951414
##  [65,]  1.2433447 -1.5554508 -9.125675e-01  0.7261706778  0.179254305
##  [66,]  0.8947850 -1.8709599 -1.091397e+00  0.6656064966  0.627036520
##  [67,]  1.0102166 -2.4878734 -1.551434e+00 -0.0136513204  0.280841967
##  [68,]  1.3020091 -2.0579575 -1.547957e+00  0.1327128719  0.166145653
##  [69,]  0.6413830 -2.3666940 -1.679636e+00  0.4737128211  0.068602792
##  [70,]  1.1486187 -2.4380403 -1.052639e+00 -0.2025685176  0.319521549
##  [71,]  0.8931504 -2.6309922 -1.133818e+00 -0.2153469785  0.306980519
##  [72,]  0.9597449 -1.9109350 -1.179982e+00  0.5290306276  0.200750574
##  [73,]  1.5930898 -1.7197211 -1.150464e+00  0.1633603334  0.102672605
##  [74,]  1.2896989 -1.7177719 -1.203659e+00  0.3210778342  0.209526987
##  [75,]  1.0601077 -1.7068550 -1.290722e+00  0.6086513959  0.157702975
##  [76,]  0.8102744 -1.8109639 -1.214487e+00  0.8102703778  0.061647569
##  [77,]  1.0189620 -1.9274434 -9.306160e-01  0.3164933985  0.176160585
##  [78,]  1.3851612 -2.5528279 -9.969352e-01 -0.4239746217  0.025645423
##  [79,]  1.2982993 -1.6552566 -1.350489e+00  0.4225025449  0.099525834
##  [80,]  1.3250985 -1.9842294 -1.395971e+00  0.2383494947 -0.001346272
##  [81,]  1.0233621 -2.6707818 -1.300869e+00 -0.0506592000 -0.197799011
##  [82,]  0.5690813 -2.0508793 -1.399208e+00  0.7572836660  0.730254901
##  [83,]  0.5869826 -2.3898224 -1.190377e+00  0.4264943703  0.114401141
##  [84,]  0.7776370 -2.2387693 -1.609509e+00  0.4355763184  0.204489690
##  [85,]  1.5364888 -1.3534210 -9.365923e-01  0.5091501588 -0.018852710
##  [86,]  0.8004368 -2.3962347 -1.292215e+00  0.1766841360  0.257357348
##  [87,]  0.9784944 -1.9494583 -1.268897e+00  0.6893943726  0.072356386
##  [88,]  2.0486946 -1.1347197 -4.613436e-01 -0.2167658268 -0.462841227
##  [89,]  2.0609374 -1.0200433 -5.568194e-01 -0.0518956951 -0.495092784
##  [90,]  1.5027492 -1.4517344 -5.914508e-01 -0.2638182862 -0.122625201
##  [91,]  1.2454190 -1.6248447 -8.034087e-01 -0.2572522574 -0.251328679
##  [92,]  1.5166008 -1.3267635 -5.471580e-01 -0.3257302533 -0.074542288
##  [93,]  1.9244223 -1.2688693 -5.310449e-01 -0.3099215179 -0.434035587
##  [94,]  1.3645560 -1.4560181 -6.251179e-01 -0.2166348120 -0.184921760
##  [95,]  1.4332766 -1.5162627 -5.820915e-01 -0.2354755698 -0.126479251
##  [96,]  0.9918570 -2.2555878 -1.613563e+00  0.2754939353  0.208568399
##  [97,]  1.5596143 -1.3902318 -6.135966e-01 -0.2792824398 -0.089462382
##  [98,]  1.4269160 -1.3723623 -5.923807e-01 -0.1516963784 -0.234637468
##  [99,]  1.4200960 -1.3633450 -5.039183e-01 -0.1969147642 -0.198652435
## [100,]  1.3137432 -1.5985442 -6.712030e-01 -0.3492518398 -0.119092744
## [101,]  2.2849363 -1.5666596 -4.046573e-01 -0.9472347549 -0.356001094
## [102,]  2.3037245 -1.4884207 -3.219346e-01 -1.0508207887 -0.266098356
## [103,]  2.3056006 -1.6145610 -3.420332e-01 -1.0498265728 -0.254167419
## [104,]  1.5016417 -1.4442199 -6.100707e-01 -0.2870215301 -0.127458659
## [105,]  1.4048626 -1.4913675 -5.833225e-01 -0.3179546343 -0.129606928
## [106,]  1.4185876 -1.5463000 -5.962801e-01 -0.2734244071 -0.123176696
## [107,]  1.4268391 -1.4478407 -4.698280e-01 -0.2641568126 -0.116341572
## [108,]  1.4837139 -1.4682829 -6.606025e-01 -0.2477492203 -0.128813668
## [109,]  1.4193955 -1.5341650 -6.486014e-01 -0.2782496839 -0.091925164
## [110,]  1.3240381 -1.7161573 -7.530511e-01 -0.3206843827 -0.273056735
## [111,]  2.1063472 -1.0464733 -4.266070e-01 -0.2704645881 -0.465663571
## [112,]  1.3342987 -1.7534853 -6.069081e-01 -0.4790528229 -0.047590734
## [113,]  2.2726790 -1.7146245 -4.887280e-01 -0.8116162403 -0.492114347
## [114,]  1.2630421 -1.6255564 -7.996530e-01 -0.1323756898 -0.353429615
## [115,]  1.8915609 -1.1622540 -6.963353e-01  0.0380762068 -0.735496810
## [116,]  2.0914790 -1.1000042 -5.711938e-01 -0.0746904832 -0.675638905
## [117,]  1.9509668 -1.1459015 -4.625413e-01 -0.0186726467 -0.599258512
## [118,]  1.9056501 -1.1391890 -4.573904e-01 -0.1081517422 -0.594096569
## [119,]  1.3851612 -2.5528279 -9.969352e-01 -0.4239746217  0.025645423
## [120,] -1.7821061  0.9357388 -3.527896e-02  0.2544845538 -0.061422205
## [121,] -1.6097056  0.6875606  7.349768e-02  0.3579917640 -0.063730262
## [122,] -1.7178189  0.8063983 -5.887644e-02  0.1790821102 -0.013018660
## [123,] -2.4837683 -0.3081262  3.145661e-02 -1.0078017509  0.414603441
## [124,] -1.8079753  0.8523505 -6.290015e-02  0.0415369064 -0.017237037
## [125,] -2.5498135 -0.3418270  4.702462e-02 -0.9636759654  0.360102595
## [126,] -1.9567671  0.7250382 -8.705428e-02 -0.0163257978  0.018757254
## [127,] -1.5287712  0.3373404  2.103437e-01  0.2777613118  0.011419799
## [128,] -1.6855350  0.7068882  6.869857e-02  0.0849181309  0.010110292
## [129,] -1.6799522  0.6535337  8.597124e-02  0.1353267481  0.003405874
## [130,] -1.7182793  0.3356592  1.020004e-01 -0.2273033831  0.080802662
## [131,] -1.8879987  0.6702338 -8.035208e-02 -0.0317601305 -0.029112603
## [132,] -1.3976740  1.0994870  4.071836e-02  0.4145454815 -0.109000799
## [133,] -1.4333008  1.2914902 -1.773336e-02  0.4965323277 -0.087440909
## [134,] -1.2593338  0.7786579  2.459925e-01  0.4450095104 -0.189308605
## [135,] -1.8299585  0.4149523  1.004233e-01 -0.0492359498  0.088868786
## [136,] -1.9316993  0.5052964  9.230805e-02 -0.2533784207  0.076782322
## [137,] -1.2593338  0.7786579  2.459925e-01  0.4450095104 -0.189308605
## [138,] -1.7250673  0.8773319 -2.039163e-02 -0.0924241557  0.096784046
## [139,] -1.7216083  0.5924819  9.725722e-02  0.0159373681  0.046611160
## [140,] -1.4335765  1.0295516 -9.708944e-03  0.2245546541 -0.162214502
## [141,] -1.4922623  0.9138837  1.316622e-01  0.2608351056 -0.026567419
## [142,] -1.3624677  1.1224589  1.963760e-01  0.5261209127 -0.087464446
## [143,] -1.3857731  1.2490095  7.209095e-02  0.3961710681 -0.160122584
## [144,] -1.4063105  1.0460328  7.614225e-02  0.3397387687 -0.211837228
## [145,] -2.2111365  0.1994839  1.514631e-01 -0.7776777846  0.217017659
## [146,] -1.7182793  0.3356592  1.020004e-01 -0.2273033831  0.080802662
## [147,] -1.8794584  0.8215182 -3.169290e-02 -0.0936000726  0.010743403
## [148,] -1.4593233  1.8490788 -2.651474e-01  0.0972620651 -0.095181423
## [149,] -1.7397189  1.4120162 -2.532970e-01 -0.1193288439  0.039640443
## [150,] -1.8208267  1.3689187 -3.266021e-01 -0.1953787987  0.022203883
## [151,] -1.2910528  1.6782206 -2.103145e-01  0.4302934905 -0.133664854
## [152,] -1.5242436  1.8044204 -3.084560e-01  0.2838971686 -0.131189851
## [153,] -1.7031635  1.3690717 -2.238618e-01 -0.2091441101  0.073266443
## [154,] -1.7782382  1.6281342 -3.379413e-01  0.0128006873 -0.067030299
## [155,] -1.9127667  1.1563632 -3.311623e-01 -0.2179915875  0.012067083
## [156,] -1.7519713  1.3991796 -2.795293e-01 -0.2310447026  0.079536758
## [157,] -1.7708516  1.4183137 -2.730061e-01  0.0607730724 -0.060473873
## [158,] -1.3379023  1.6308089 -2.113134e-01  0.2506588925 -0.081752927
## [159,] -1.8259271  1.4764561 -2.970204e-01 -0.0862932946  0.003479588
## [160,] -1.2634364  1.8231629 -1.440486e-01  0.1718464207  0.138925797
## [161,] -1.7706674  1.4260491 -2.720666e-01 -0.1182539818  0.027514737
## [162,] -1.7550485  1.8266947 -4.228753e-01  0.1966691954 -0.084424390
## [163,] -1.5740527  1.7736688 -3.197549e-01  0.1655929520 -0.074360762
## [164,] -1.3334085  1.7717299 -2.737830e-01  0.3743999106 -0.127775538
## [165,] -1.8621680  1.4065277 -3.433491e-01 -0.1690710738  0.026605966
## [166,] -1.4817457  2.0361137 -3.795028e-01  0.1908833912 -0.077017356
## [167,] -1.8240038  1.2925430 -2.615833e-01 -0.2066553255  0.037351493
## [168,] -1.9036350  1.1927811 -2.861868e-01 -0.3530678431  0.038456663
## [169,] -2.0020922  1.2314371 -3.184723e-01 -0.4074384473  0.080822078
## [170,] -1.6863456  1.5552550 -3.163020e-01  0.0260960867 -0.026001557
## [171,] -1.6857115  1.4701915 -2.547130e-01 -0.1351827805 -0.119117340
## [172,] -1.5645926  1.6693219 -2.820254e-01  0.1090094869 -0.206311546
## [173,] -1.6431149  1.6369257 -2.698690e-01 -0.0003442544 -0.049645484
## [174,] -1.6165737  1.6707067 -2.985345e-01  0.0885303224 -0.037445380
## [175,] -1.5160929  1.8100963 -2.736298e-01  0.1831343281 -0.105466106
## [176,] -1.4291019  1.8557958 -2.503902e-01  0.2488641895 -0.107735354
## [177,] -1.7623546  1.2124241 -2.232618e-01 -0.3032582183 -0.170488901
## [178,] -1.7770965  1.6334597 -4.156476e-01 -0.0760106644  0.012877149
## [179,] -1.8539682  1.5772551 -4.342079e-01 -0.0950046531 -0.002425308
## [180,] -1.9563477  1.2056229 -3.110391e-01 -0.1550888312  0.013671978
## [181,] -1.8339199  1.5913197 -4.529597e-01 -0.0907722155 -0.004327614
## [182,] -1.3812960  1.7937448 -3.439181e-01  0.4492210801 -0.126572656
## [183,] -1.7339470  1.6976485 -3.882368e-01  0.1142546685 -0.009199461
## [184,] -1.7058367  1.5993365 -3.848950e-01 -0.0014980557  0.004726687
## [185,] -1.6597247  1.7015571 -3.849047e-01 -0.1817934288  0.027704203
## [186,] -1.7509256  1.6772132 -4.586073e-01  0.0784618988 -0.053039296
## [187,] -1.6946457  1.6702985 -4.181720e-01 -0.0052462430 -0.004766333
## [188,] -1.8690746  1.5451549 -4.645338e-01 -0.2354901517  0.027750525
## [189,] -1.7709664  1.4726949 -5.213314e-01 -0.3099356435 -0.165587633
## [190,] -1.7897240  1.5398215 -4.397458e-01 -0.0656828124 -0.010422426
## [191,] -2.0742731  0.8833671 -1.722746e-01 -0.2645846386  0.110433276
## [192,] -1.3904139  2.1076960 -4.205864e-01  0.3108469297 -0.135178806
## [193,] -1.7377423  1.6980129 -4.950300e-01 -0.1707365346 -0.002580614
## [194,] -1.8323276  1.6578070 -4.514545e-01 -0.1979256139  0.026239523
## [195,] -1.8672861  1.6777345 -4.556293e-01 -0.0399614874 -0.014226888
## [196,] -1.6088333  1.4289516 -2.855827e-01  0.2124194163 -0.059425946
## [197,] -1.8241903  1.4365231 -3.900331e-01 -0.2318969126  0.067676333
## [198,] -1.5488523  1.6818865 -3.437662e-01  0.0581391259 -0.033664143
## [199,] -1.7022993  1.6636683 -4.876881e-01 -0.3490398234 -0.154313366
## [200,] -1.9094050  1.5220593 -4.046353e-01 -0.2913380635  0.067649162
## [201,] -1.6406553  1.9871102 -4.740963e-01  0.0825407704 -0.087342860
## [202,] -1.7071054  1.9787327 -4.996512e-01  0.2436043361 -0.121964525
## [203,] -1.7908242  1.6611027 -4.463550e-01 -0.1363032910 -0.007900480
## [204,] -1.6289080  2.0826485 -5.711895e-01 -0.0460860241 -0.015349914
## [205,] -1.9903369  1.3222111 -4.307109e-01 -0.4140262008  0.077044732
## [206,] -1.7267181  1.6174648 -3.560766e-01 -0.0211010840 -0.022537063
## [207,] -2.0511518  0.4407758 -1.239556e-01 -0.1654892867  0.054734780
## [208,] -1.4187340  0.9955696 -7.028759e-03  0.5097265693 -0.086876308
## [209,] -1.7244772  0.7989221  5.901215e-02  0.0872590764  0.061267151
## [210,] -2.4572491  0.1595509 -2.868641e-01 -0.7815829395  0.243614030
## [211,] -1.6751303  0.7614038 -4.316831e-02  0.2262013713 -0.050605850
## [212,] -1.7126692  0.7752786 -3.723766e-02  0.3069598594 -0.017574873
## [213,] -2.5816073 -0.4652796  6.309360e-02 -0.7813544450  0.323301541
## [214,] -1.6279961  0.7537831 -1.844262e-02  0.1306665917 -0.350779579
## [215,] -1.8115248  1.1063091 -7.734399e-02  0.1190940094 -0.031264854
## [216,] -1.8554293  0.9253740 -1.588235e-01  0.1042746592  0.021987080
## [217,] -2.6481066 -0.2787996 -7.235941e-02 -0.8722492036  0.343955004
## [218,] -1.5483463  1.1834618 -5.198809e-02  0.3982239133 -0.149425102
## [219,] -1.7630072  0.5600795 -1.055403e-01  0.1994110174  0.017985654
## [220,] -1.9151299  0.6613004 -1.338789e-02  0.0717021199  0.109138646
## [221,] -1.9732553  0.4364541 -6.488793e-02 -0.0815661887  0.106986872
## [222,] -1.9271485  0.6438785 -6.922989e-02 -0.0325491686  0.088619098
## [223,] -1.8582907  0.5475690  3.258020e-02  0.1356881583  0.029158570
## [224,] -1.7459325  0.7811679  3.170352e-02  0.1304458635 -0.043836837
## [225,] -1.7743513  0.9874912 -7.619610e-02  0.1706258319  0.018918311
## [226,] -1.8564759  0.6412486 -1.210501e-02 -0.0118052550  0.010652649
## [227,] -1.8215076  0.6778087  7.487880e-03 -0.0018364731  0.055289102
## [228,] -1.4726670  0.9459139  8.328993e-02  0.5369939397 -0.124044993
## [229,] -1.5613641  0.7512569  2.992360e-02  0.2238667572  0.037671605
## [230,] -1.8342202  0.6924831  6.625196e-03 -0.0768086714  0.094293955
## [231,] -1.5384088  0.7392862  7.866684e-02  0.0590355141 -0.223418689
## [232,] -1.7328789  0.7193759  5.147681e-02  0.1530394693  0.014006260
## [233,] -1.6291117  0.7972391  1.202807e-01  0.5505075044 -0.023930604
## [234,] -1.7573656  0.6446385  8.172953e-03 -0.0352054622 -0.184397063
## [235,] -1.9891863  0.5021243 -2.138313e-02 -0.1201205572  0.124730551
## [236,] -1.6477838  1.0438009  3.632101e-02  0.4056517096 -0.092016115
## [237,] -2.1823786  0.2621524 -2.198562e-02 -0.5691155558  0.057489251
## [238,] -1.8863207  0.5862007 -3.238206e-02 -0.0683081684  0.024090434
## [239,] -1.8215076  0.6778087  7.487880e-03 -0.0018364731  0.055289102
## [240,] -1.0562819 -1.8849586  8.134634e-01 -0.1009130167  0.080734261
## [241,] -1.1718574 -1.8078398  8.061754e-01  0.0987556257  0.022996995
## [242,] -1.1747986 -1.7751811  7.308878e-01  0.0160277466  0.054606827
## [243,] -1.0788743 -1.8125005  7.001711e-01  0.1734592492  0.047930806
## [244,] -1.0790264 -1.5817939  7.187831e-01  0.1276327602  0.057276047
## [245,] -1.2694406 -2.0712740  6.550343e-01 -0.0916768440  0.154884201
## [246,] -1.4343620 -2.0603104  9.726186e-01 -0.0427461898  0.111587294
## [247,] -1.1775638 -2.0507147  8.706037e-01  0.1158017679  0.027574258
## [248,] -1.2277190 -1.8525166  8.516937e-01 -0.1421591246  0.044765119
## [249,] -0.7622691 -1.6803033  7.506781e-01  0.0748230386  0.066874442
## [250,] -1.0602926 -1.8637091  7.761526e-01  0.0696483051  0.007015509
## [251,] -1.0045531 -1.7306684  5.832766e-01  0.0586719111  0.014901884
## [252,] -1.0657651 -2.0106009  7.477518e-01 -0.0866299169  0.100240658
## [253,] -1.0337915 -1.8980999  7.589126e-01  0.0779803537  0.108164405
## [254,] -1.1251976 -2.0228537  8.143071e-01 -0.0912732448  0.093620108
## [255,] -1.1469355 -1.9025984  9.034707e-01  0.1614667314  0.054480627
## [256,] -1.1678685 -1.8356307  8.559839e-01 -0.0739556183  0.132748873
## [257,] -1.2113990 -1.9055903  7.778266e-01 -0.0039203987  0.074493345
## [258,] -1.1994321 -1.8873848  6.965410e-01 -0.0696616621  0.117271583
## [259,] -1.1234870 -1.6904939  6.192902e-01 -0.0426787868  0.074839230
## [260,] -1.0009875 -1.9358364  6.180061e-01 -0.0192021957  0.041898530
## [261,] -1.0653296 -1.8722276  6.888437e-01  0.0966835108  0.010058996
## [262,] -1.1340498 -1.8698826  7.378555e-01  0.0521185808  0.033490320
## [263,] -1.1264476 -1.8155883  6.761190e-01  0.1115936289  0.055738624
## [264,] -1.3906598 -2.2635686  9.332784e-01 -0.0956292166  0.093748967
## [265,] -1.0700924 -1.8816391  7.735862e-01  0.0535670529  0.099686512
## [266,] -1.1850105 -1.8157484  7.388524e-01  0.0625574291  0.014437698
## [267,] -1.1291235 -1.8039368  6.518559e-01 -0.0471697274  0.044120573
## [268,] -1.1607513 -1.9235745  8.877838e-01 -0.0136013556  0.088103463
## [269,] -0.6162255 -0.6509414  5.701714e-01 -0.2669313379  0.037726747
## [270,] -0.5943386 -0.3976471  3.524611e-01 -0.2484317278  0.025685452
## [271,] -0.8636353 -0.9202528  4.929314e-01 -0.1271878424 -0.024301448
## [272,] -0.8943741 -0.7665976  5.292233e-01 -0.1140331007 -0.026369333
## [273,] -0.9406739 -0.6569263  4.539243e-01 -0.0018501865 -0.068884633
## [274,] -0.9831377 -0.9716010  5.039806e-01 -0.1422556128 -0.029824184
## [275,] -0.6422548 -0.5815489  4.805844e-01 -0.0442746257 -0.002984974
## [276,] -0.6816623 -0.6933229  4.888584e-01 -0.0217855876  0.073367199
## [277,] -0.5963722 -0.6275245  4.743971e-01 -0.0330329530 -0.007735678
## [278,] -0.7851692 -0.9221646  6.837651e-01 -0.0729860579 -0.005213343
## [279,] -0.5995188 -0.6837602  1.253924e-01 -0.1769133862  0.170425062
## [280,] -0.8438760 -0.7946118  4.433999e-01 -0.3908144821  0.119973888
## [281,] -0.6720668 -0.8080937  5.632240e-01 -0.0760354377  0.008486348
## [282,] -0.6049674 -0.7435461  4.163869e-01 -0.2179192006  0.033789134
## [283,] -0.6750422 -0.5776411  3.998099e-01 -0.1628786172 -0.052170982
## [284,] -0.7694110 -0.8165060  5.570037e-01 -0.2467109699  0.094392357
## [285,] -0.2760655 -0.1671596  4.331238e-01 -0.0215426784 -0.040089519
## [286,] -0.5727469 -0.8238554  6.690880e-01 -0.2256065927  0.020008405
## [287,] -0.6444701 -0.4916469  4.545389e-01 -0.0555261120 -0.020104884
## [288,] -0.6449905 -0.4931173  4.619264e-01 -0.2711204336  0.039537296
## [289,] -0.6301577 -0.8957997  5.814457e-01 -0.1624954468  0.064121728
## [290,] -0.7547400 -0.5820306  4.691199e-01 -0.0719888815 -0.022974937
## [291,] -0.7132627 -0.6449542  5.064120e-01 -0.1201080302 -0.006069777
## [292,] -0.7749361 -0.8410832  5.707717e-01 -0.1913850493  0.047588587
## [293,] -0.9646953 -0.4898548  4.903633e-01 -0.0098541434 -0.051260454
## [294,] -0.6440462 -0.5869810  5.502600e-01  0.0479134289 -0.073813852
## [295,] -0.3648100 -0.3356016  4.075140e-01 -0.1124916820 -0.057131928
## [296,] -0.5346165 -0.5299572  4.257782e-01 -0.2288754481  0.030964030
## [297,] -0.3390696 -0.2428245  2.807053e-01 -0.0641818802  0.069547775
## [298,] -0.6453544 -0.5145738  3.697604e-01 -0.2488782344  0.043484390
## [299,] -0.8636353 -0.9202528  4.929314e-01 -0.1271878424 -0.024301448
## [300,] -0.8943741 -0.7665976  5.292233e-01 -0.1140331007 -0.026369333
##                  PC6           PC7
##   [1,] -2.334931e-03  9.586381e-04
##   [2,]  2.946986e-03  1.020335e-03
##   [3,]  5.414628e-03  8.303696e-04
##   [4,]  5.617206e-03  9.448638e-04
##   [5,]  1.685765e-03  9.668578e-04
##   [6,]  5.184239e-03  9.775808e-04
##   [7,]  2.062019e-03  9.137155e-04
##   [8,] -8.239978e-04  9.452848e-04
##   [9,] -2.899542e-03  8.552953e-04
##  [10,]  2.239517e-03  8.358339e-04
##  [11,]  1.245397e-02  1.050927e-03
##  [12,] -1.275477e-03  7.918359e-04
##  [13,]  4.554520e-03  7.631647e-04
##  [14,] -2.896384e-02 -2.669892e-02
##  [15,]  8.222939e-04  8.658250e-04
##  [16,]  4.738489e-03  9.420006e-04
##  [17,] -1.852536e-03  8.634677e-04
##  [18,] -1.111859e-03  7.367465e-04
##  [19,]  4.595972e-03  7.816179e-04
##  [20,]  4.528364e-03  9.173402e-04
##  [21,]  5.191144e-03  8.558793e-04
##  [22,]  2.176916e-03  9.157980e-04
##  [23,]  1.087007e-03  8.746648e-04
##  [24,]  8.161999e-04  7.297943e-04
##  [25,] -2.097494e-03  8.007140e-04
##  [26,]  3.510400e-03  8.825632e-04
##  [27,]  4.930453e-03  9.092009e-04
##  [28,]  9.360605e-04  8.064522e-04
##  [29,]  3.837552e-03  9.786595e-04
##  [30,] -6.152930e-03  4.444485e-04
##  [31,] -3.443607e-03  4.205435e-04
##  [32,] -4.771035e-03  2.762183e-04
##  [33,] -1.443949e-03  4.012768e-04
##  [34,] -4.520045e-04  3.887818e-04
##  [35,] -1.482102e-03  4.667390e-04
##  [36,]  6.258954e-04  4.333274e-04
##  [37,] -4.260919e-03  1.823387e-04
##  [38,] -6.385345e-03  2.021854e-04
##  [39,] -7.312478e-03  3.779628e-04
##  [40,]  2.585708e-03  4.460757e-04
##  [41,] -5.513992e-03  3.250956e-04
##  [42,] -7.711213e-03  4.497519e-04
##  [43,]  1.500017e-04  4.861367e-04
##  [44,]  1.503757e-03  5.253248e-04
##  [45,] -3.727092e-03  3.773450e-04
##  [46,] -4.547758e-03  2.793730e-04
##  [47,]  4.833608e-03  4.197832e-04
##  [48,]  3.446416e-03  5.521904e-04
##  [49,] -4.811751e-03  6.038312e-04
##  [50,]  4.591192e-03  4.764630e-04
##  [51,] -8.986013e-03  3.242563e-04
##  [52,] -1.898150e-03  4.294056e-04
##  [53,] -3.424880e-03  4.952911e-04
##  [54,] -4.703082e-03  4.406668e-04
##  [55,]  3.218648e-03  5.111696e-04
##  [56,] -1.472706e-03  3.963982e-04
##  [57,]  2.866702e-03  5.503168e-04
##  [58,] -3.420656e-03  4.399337e-04
##  [59,] -6.413719e-04  2.819940e-04
##  [60,] -3.100733e-03  3.093796e-04
##  [61,] -2.573234e-03  7.370832e-05
##  [62,]  6.397883e-04  1.119815e-04
##  [63,]  2.173916e-03  1.578885e-04
##  [64,] -4.181896e-03  1.282820e-04
##  [65,]  7.568126e-03  3.268089e-04
##  [66,] -9.049039e-05  4.432851e-04
##  [67,]  4.139167e-03  2.455538e-04
##  [68,] -4.404285e-03 -4.248726e-03
##  [69,]  5.003361e-03  1.993531e-05
##  [70,]  1.649970e-03  3.389591e-04
##  [71,] -3.972698e-03  1.926221e-04
##  [72,] -2.375043e-03  1.094616e-04
##  [73,]  3.518857e-03  2.427891e-04
##  [74,]  4.799009e-03  2.866220e-04
##  [75,] -8.840448e-04  9.712172e-05
##  [76,] -1.556887e-03 -1.332355e-05
##  [77,] -2.295512e-03  1.539510e-04
##  [78,] -4.081178e-03  8.538354e-05
##  [79,] -3.021695e-03  5.189815e-05
##  [80,]  7.037333e-03  1.493088e-04
##  [81,] -3.060211e-03 -1.701159e-04
##  [82,] -4.661064e-04  4.149171e-04
##  [83,]  2.432497e-03  8.266767e-05
##  [84,] -4.448334e-04  5.586236e-05
##  [85,]  5.810299e-03  2.077777e-04
##  [86,] -4.075116e-03  9.917964e-05
##  [87,] -3.462200e-03 -2.291334e-05
##  [88,] -9.274447e-04 -2.695051e-05
##  [89,] -3.116284e-03 -1.085219e-04
##  [90,] -3.939881e-04  1.361835e-04
##  [91,]  4.828714e-03  6.756417e-05
##  [92,] -5.909315e-04  1.826814e-04
##  [93,]  2.055909e-03  2.193078e-05
##  [94,]  2.046863e-03  1.097123e-04
##  [95,]  2.289092e-03  1.684205e-04
##  [96,] -2.846058e-03  5.220853e-05
##  [97,] -5.902476e-03  7.306429e-05
##  [98,]  2.573437e-03  9.296179e-05
##  [99,] -5.571797e-03 -7.546973e-07
## [100,]  3.523064e-04  1.208292e-04
## [101,] -4.940791e-03  5.691995e-05
## [102,]  1.490353e-04  2.289201e-04
## [103,] -6.752324e-03  1.161451e-04
## [104,] -3.931862e-04  1.316796e-04
## [105,]  3.200727e-03  1.848394e-04
## [106,]  5.574712e-03  2.235759e-04
## [107,]  1.738856e-03  1.870901e-04
## [108,] -2.108492e-03  8.889964e-05
## [109,] -3.982094e-03  7.902103e-05
## [110,] -3.391575e-03 -6.586225e-05
## [111,] -4.978449e-03 -7.850592e-05
## [112,] -3.780373e-03  1.196860e-04
## [113,] -3.511529e-04  9.969070e-06
## [114,]  1.966761e-03 -5.802430e-05
## [115,] -2.169249e-04 -2.783924e-04
## [116,]  2.002938e-03 -1.489240e-04
## [117,]  4.062906e-03 -6.458508e-05
## [118,]  2.693853e-03 -8.108225e-05
## [119,] -4.081178e-03  8.538354e-05
## [120,] -8.993698e-04 -6.815043e-05
## [121,] -2.672885e-03 -7.791517e-05
## [122,]  2.135213e-04 -1.054438e-05
## [123,] -1.553842e-03  2.485740e-04
## [124,] -1.636843e-03 -4.328676e-05
## [125,] -7.553487e-03  1.019936e-04
## [126,]  4.142685e-03  5.923808e-05
## [127,]  1.580272e-03  7.203382e-05
## [128,]  8.019282e-03  1.630935e-04
## [129,]  7.824975e-03  1.534025e-04
## [130,] -4.039928e-03  2.664346e-05
## [131,]  4.844142e-03  4.553493e-05
## [132,]  3.657735e-04 -3.485138e-05
## [133,] -9.417009e-04 -5.542300e-05
## [134,]  6.677950e-03  5.025862e-05
## [135,] -6.385729e-04  6.579714e-05
## [136,]  4.210125e-03  1.431292e-04
## [137,]  6.677950e-03  5.025862e-05
## [138,] -2.396517e-03  4.963667e-05
## [139,]  4.808178e-03  1.378310e-04
## [140,]  7.447497e-03  4.607483e-05
## [141,]  7.882814e-03  1.582472e-04
## [142,]  2.064825e-03  2.900754e-05
## [143,]  3.714569e-03 -3.438925e-06
## [144,] -9.396063e-04 -1.197403e-04
## [145,]  4.548882e-04  1.882905e-04
## [146,] -4.039928e-03  2.664346e-05
## [147,]  1.423652e-03  3.346475e-05
## [148,]  9.830711e-04 -2.656776e-05
## [149,] -1.562628e-03  1.575088e-06
## [150,]  1.042417e-02  1.730880e-04
## [151,]  3.366550e-03 -1.603440e-05
## [152,]  4.438447e-03 -2.271053e-05
## [153,]  1.084609e-03  8.302135e-05
## [154,]  7.786497e-03  6.125300e-05
## [155,]  6.341227e-03  8.365973e-05
## [156,]  7.807913e-03  1.877148e-04
## [157,]  9.401088e-03  9.463983e-05
## [158,]  9.333272e-03  1.261805e-04
## [159,] -1.237076e-03 -3.502092e-05
## [160,]  9.154033e-03  3.066016e-04
## [161,]  3.657631e-03  7.399152e-05
## [162,] -9.834140e-03 -1.188457e-02
## [163,]  2.629399e-03 -1.236017e-05
## [164,]  5.546849e-03  1.643078e-05
## [165,] -2.914883e-02  2.662694e-02
## [166,] -2.988035e-01 -5.022662e-03
## [167,]  1.336811e-04  2.114655e-05
## [168,] -1.761130e-03 -1.387666e-05
## [169,] -2.921433e-03 -1.415106e-05
## [170,]  1.276226e-03 -8.068291e-06
## [171,]  6.676515e-03  3.592754e-05
## [172,] -3.762652e-04 -1.462409e-04
## [173,] -1.869376e-03 -6.097852e-05
## [174,]  4.209770e-03  4.142128e-05
## [175,]  9.115410e-03  8.686062e-05
## [176,]  7.159470e-03  6.205511e-05
## [177,]  4.131011e-03 -4.005715e-05
## [178,]  3.907947e-03  4.715226e-05
## [179,] -2.061634e-03 -7.402033e-05
## [180,]  6.621688e-03  8.464812e-05
## [181,] -1.771250e-03 -7.125222e-05
## [182,]  1.775980e-03 -6.639778e-05
## [183,]  6.615772e-03  7.346921e-05
## [184,]  3.755776e-03  4.500923e-05
## [185,]  5.154454e-03  1.050580e-04
## [186,]  9.080078e-03  7.309759e-05
## [187,]  1.102714e-03 -7.630574e-06
## [188,]  4.105769e-03  5.269659e-05
## [189,]  3.725076e-03 -8.359096e-05
## [190,]  3.345786e-03  1.320786e-05
## [191,] -3.762297e-03 -1.181886e-05
## [192,]  3.836173e-03 -3.279601e-05
## [193,]  7.840782e-04 -1.522551e-05
## [194,]  6.037409e-03  8.994462e-05
## [195,]  8.973716e-03  9.514429e-05
## [196,] -8.857651e-04 -7.119852e-05
## [197,]  8.660213e-03  1.696731e-04
## [198,]  8.525194e-03  1.185052e-04
## [199,] -9.166599e-04 -1.324960e-04
## [200,]  3.997787e-03  8.730080e-05
## [201,]  7.303820e-03  3.673733e-05
## [202,]  3.592087e-03 -7.222833e-05
## [203,]  4.729065e-03  4.514915e-05
## [204,] -5.966879e-04 -4.642964e-05
## [205,] -1.473814e-03 -6.143030e-06
## [206,]  1.773602e-03 -2.828798e-06
## [207,]  5.708078e-03  9.801102e-05
## [208,] -1.844410e-03 -7.552275e-05
## [209,] -2.947616e-03  1.298256e-05
## [210,] -5.181870e-03  1.759866e-05
## [211,] -8.372866e-04 -2.204029e-03
## [212,]  1.718573e-03  5.104281e-06
## [213,]  1.302333e-03  2.061922e-04
## [214,]  6.984948e-04 -2.203949e-04
## [215,]  7.158750e-03  9.157450e-05
## [216,]  7.754585e-03  1.176175e-04
## [217,] -2.471434e-03  1.410919e-04
## [218,]  5.610749e-03 -3.019019e-06
## [219,]  4.196640e-03  5.733116e-05
## [220,]  3.206380e-03  1.146765e-04
## [221,]  1.557680e-03  7.678795e-05
## [222,]  7.712656e-04  4.873970e-04
## [223,]  3.918714e-03  7.625577e-05
## [224,]  8.178222e-03  1.144514e-04
## [225,]  7.705444e-03  1.330061e-04
## [226,]  7.785653e-03  1.340595e-04
## [227,] -3.558864e-04  3.684021e-05
## [228,]  8.265827e-03  7.175954e-05
## [229,]  7.761227e-03  1.762725e-04
## [230,] -3.866410e-03  1.005393e-05
## [231,]  7.607348e-03  1.307489e-05
## [232,]  8.094905e-03  1.546949e-04
## [233,]  5.748314e-03  8.404796e-05
## [234,]  8.176009e-03  1.981917e-05
## [235,]  1.780728e-03  1.025228e-04
## [236,]  1.349763e-03 -3.478689e-05
## [237,]  5.871209e-03  1.293331e-04
## [238,]  1.305622e-03  3.198549e-05
## [239,] -3.558864e-04  3.684021e-05
## [240,] -1.932413e-03  1.761227e-04
## [241,] -7.527321e-03 -4.137051e-04
## [242,] -6.460027e-03  5.136341e-05
## [243,] -3.888390e-03  8.245441e-05
## [244,] -6.393284e-04  1.552970e-04
## [245,]  9.210361e-04  2.223481e-04
## [246,] -5.533328e-03  1.135069e-04
## [247,] -3.325396e-03  9.167809e-05
## [248,] -7.877363e-03  4.365390e-05
## [249,]  2.728503e-03  2.580585e-04
## [250,]  3.850058e-03  2.026722e-04
## [251,] -3.903406e-03  5.929364e-05
## [252,] -8.097221e-04  1.928150e-04
## [253,]  3.132178e-03  2.600207e-04
## [254,]  2.584233e-03  2.487256e-04
## [255,] -9.737166e-04  1.585437e-04
## [256,] -2.815556e-03  1.918429e-04
## [257,] -4.552084e-03  9.846701e-05
## [258,] -3.484038e-03  1.399584e-04
## [259,]  3.081216e-03  2.187307e-04
## [260,] -1.394134e-03  1.260076e-04
## [261,]  4.161266e-02  8.165585e-04
## [262,] -7.268287e-03  2.349265e-05
## [263,]  1.252997e-03  1.690255e-04
## [264,] -6.471700e-03  8.248337e-05
## [265,] -7.051018e-03  8.534456e-05
## [266,]  1.622405e-03  1.534450e-04
## [267,] -2.163787e-03  1.118729e-04
## [268,] -1.671471e-03  1.787121e-04
## [269,]  1.006368e-03  2.482562e-04
## [270,]  4.400769e-03  2.698942e-04
## [271,]  1.819712e-03  1.627628e-04
## [272,]  1.015797e-04  1.381950e-04
## [273,] -3.605187e-03  2.482330e-05
## [274,]  3.382974e-03  1.736446e-04
## [275,] -5.053270e-03  8.783070e-05
## [276,]  4.592830e-03  2.941501e-04
## [277,]  2.518540e-03  2.123635e-04
## [278,] -6.498058e-03  7.193468e-05
## [279,]  4.142761e-03  3.179813e-04
## [280,] -2.524562e-03  2.081619e-04
## [281,]  3.683922e-03  2.473004e-04
## [282,] -2.106245e-03  1.649015e-04
## [283,]  3.707313e-03  1.918155e-04
## [284,]  1.478476e-03  2.714968e-04
## [285,] -4.213680e-03  1.166953e-04
## [286,] -6.107405e-03  1.301084e-04
## [287,] -5.743411e-04  1.492376e-04
## [288,] -3.547164e-03  1.582288e-04
## [289,] -4.016922e-03  1.693975e-04
## [290,] -5.606027e-03  5.289041e-05
## [291,] -2.938718e-04  1.651740e-04
## [292,]  4.586480e-03  2.873838e-04
## [293,] -2.839315e-03  5.796378e-05
## [294,]  1.068278e-04  1.280177e-04
## [295,]  5.072415e-03  2.480063e-04
## [296,]  5.073974e-02 -3.466824e-02
## [297,] -8.207324e-05  2.326022e-04
## [298,]  2.106968e-03  2.384104e-04
## [299,]  1.819712e-03  1.627628e-04
## [300,]  1.015797e-04  1.381950e-04

Так мы вытаскиваем сами значения главных компонент для всех наблюдений. Их далее можем использовать, например, в регрессии.

А что еще можно сделать с самими результатами PCA? Можем посмотреть, как соотносятся категориальные переменные с главными компонентами.

pizza_pca <- tibble(PC1 = pca$x[,1], PC2 = pca$x[,2], brand = pizza$brand)
ggplot(pizza_pca, aes(PC1, PC2, color = brand)) +
  geom_point()

Вау! У нас тут даже явные кластеры нарисовались! А вот как бы это интерпретировать?

Здесь нам понадобиться матрица нагрузок (собственно, корреляционная матрица, которую мы видели в аутпуте PCA):

pca$rotation[,1:2] # сразу вытащим только первые две компоненты

##                PC1        PC2
## mois    0.06470937 -0.6282759
## prot    0.37876090 -0.2697067
## fat     0.44666592  0.2343791
## ash     0.47188953 -0.1109904
## sodium  0.43570289  0.2016617
## carb   -0.42491371  0.3203121
## cal     0.24448730  0.5674576

И поставим рядом описание переменных, чтобы было удобнее интерпретировать:

• brand — Pizza brand (class label)
• id — Sample analysed
• mois — Amount of water per 100 grams in the sample
• prot — Amount of protein per 100 grams in the sample
• fat — Amount of fat per 100 grams in the sample
• ash — Amount of ash per 100 grams in the sample
• sodium — Amount of sodium per 100 grams in the sample
• carb — Amount of carbohydrates per 100 grams in the sample
• cal — Amount of calories per 100 grams in the sample

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

Один из вариантов интерпретации: PC1 — пищевая ценность, PC2 — энергетическая ценность. Можно придумать и другие, но хорошо бы чтобы они были осмысленные.

Ну, а теперь, понимая содержания компонент, можно и на предыдущую диаграмму рассеяния осмысленно посмотреть.

Эксплораторный факторный анализ (exploratory factor analysis)

Математическая модель эскплораторного факторного анализа

У нас есть все та же матрица наблюдений \(\boldsymbol{X}^\mathrm{T}= \begin{pmatrix} \boldsymbol{X}_1 & \boldsymbol{X}_2 & \dots & \boldsymbol{X}_k \end{pmatrix}\), только мы ее транспонировали для будущего удобства. Мы предполагаем, что под нашими наблюдениями спрятаны некие факторы в количестве \(p\) штук, \(p < k\). Их мы также можем оформить в матрицу \(\boldsymbol{F}^\mathrm{T}= \begin{pmatrix} \boldsymbol{F}_1 & \boldsymbol{F}_2 & \dots & \boldsymbol{F}_p \end{pmatrix}\). Эти факторы объясняют имеющиеся переменные. Делают это они следующим образом:

\[ \boldsymbol{X}_i = a_{i1} \boldsymbol{F}_1 + a_{i2} \boldsymbol{F}_2 + \dots + a_{ip} \boldsymbol{F}_p + \boldsymbol{U}_i, \; 1 \leq i \leq k \]

\[ \boldsymbol{X}= \boldsymbol{A}\boldsymbol{F}+ \boldsymbol{U}, \\ \boldsymbol{A}= (a_{ij}), 1 \leq i \leq k, 1 \leq j \leq p, \\ \boldsymbol{U}^\mathrm{T}= \begin{pmatrix} \boldsymbol{U}_1 & \boldsymbol{U}_2 & \dots & \boldsymbol{U}_k \end{pmatrix} \]

Здесь \(\boldsymbol{U}\) — то, что не удалось объяснить факторами (остатки, уникальность, uniqueness).

Также выдвигаем ряд дополнительных предположения для формализации требования и упрощения вычислений:

1. \(\mathbb{E}\boldsymbol{X}= 0\). Так сделать можно, так как мы все равно будем стандартизировать переменные, и математическое ожидание обратиться в ноль.
2. \(\mathrm{corr}(\boldsymbol{F}_j, \boldsymbol{F}_t) = 0, \forall j \; \forall t, j \neq t, 1 \leq j \leq p, 1 \leq t \leq p\). То есть хотим, чтобы факторы были независимы (некоррелированы).
3. \(\mathrm{var}(\boldsymbol{F}) = \boldsymbol{I}\). Если «истинная» дисперсия факторов будет отличаться от единицы, то разница уйдет в матрицу \(\boldsymbol{A}\).
4. \(\mathrm{corr}(\boldsymbol{U}_i, \boldsymbol{U}_r) = 0; \; \mathrm{corr}(\boldsymbol{U}_i, \boldsymbol{F}_j) = 0, \forall i \; \forall r\; \forall j, i \neq r, 1 \leq i \leq k, 1 \leq r \leq k ,1 \leq j \leq p\). Ну, они же все-таки уникальности.

Элементы матрицы \(\boldsymbol{A}\) называются факторными нагрузками (factor loadings).

Элементы вектора \(\boldsymbol{U}\) называются уникальными факторами (specific variates).

Теперь об информативности. Аналогично PCA, возьмем в качестве меры информативности дисперсию. Мы хотим узнать, насколько хорошо факторы объясняют исходные переменные. Дисперсия переменных будет складываться из следующего:

\[ \mathrm{var}(\boldsymbol{X}_i) = \sum_{j=1}^p a_{ij}^2 + \mathrm{var}(\boldsymbol{U}_i) \]

То есть, чем больше уникальность (uniqueness) — часть дисперсии переменной, объясненной уникальными факторами — тем хуже наши факторы объясняют переменную. Что делать? Либо подобрать другую модель (изменить количество факторов), либо не брать такую переменную в факторный анализ…

Однако в факторном анализе есть еще одна интересная деталь. Мы работаем с выборочными корреляциями. По этой причине мы не можем подобрать единственно верное решение задачи факторного анализа. Матрица факторных нагрузок определена до ортогонального преобразования (по-русски, вращения.) То есть у нас есть множество решений, которые отличаются друг от друга поворотом, и мы можем выбирать самое «симпатишное».

За такую неопределенность факторный анализа часто критикуют. Да и вообще, как вы могли заметить, в факторном анализе субъективность чуть ли не на каждом шагу.

Вращение факторов (factor rotation)

Вращать вообще-то не обязательно. Если хочется, чтобы все было похоже на PCA, то можно ввести дополнительное условие, чтобы матрица \(\boldsymbol{A}^\mathrm{T}\mathrm{diag}(\mathrm{var}(\boldsymbol{U})) \boldsymbol{A}\) быда диагональной и элементы на её главной диагонали стояли в убывающем порядке.

Но, господи, вращать это же весело!

Поэтому придумали много всяких разных вариантов вращения. Самый популярный — varimax. Его идея состоит в том, чтобы найти наиболее «контрастное» решений — чтобы значения факторов были или большими, или маленькими. Эта логика хорошо согласуется с тем, что мы хотели, чтобы каждый фактор описывал что-то своё. Попутно мы минимизируем количество переменных, которые имеют высокую нагрузку на каждый фактор. Работает varimax по слудеюшей математической модели:

\[ \begin{cases} a^2 + b^2 \rightarrow \max \\ a + b = 1 \end{cases} \]

Коротенечко о других вращениях:

• quartimax — минимизирует количество факторов, необходимых для объяснения каждой переменной
• equamax — количество переменных, сильно нагружающих фактор, и количество факторов, необходимых для объяснения переменных, минимальны
• promax — наклонное вращение, позволяющее коррелировать факторы

Реализация EFA в R

Работаем с тем же датасетом. Для факторного анализа есть функция factanal(). Она ждет от нас данные и количество факторов.

factanal(pizza2, factors = 5)

## Error in factanal(pizza2, factors = 5): 5 factors are too many for 7 variables

Опа, попросили слишком много факторов. Действительно, если математические ограничения на количество факторов, оно завязано на число имеющихся переменных. Но нам не надо об этом париться, R за нас сам все проверит и посчитает.

Попробуем 4:

factanal(pizza2, factors = 4)

## Error in factanal(pizza2, factors = 4): 4 factors are too many for 7 variables

Тоже много. Давайте 3:

fan <- factanal(pizza2, factors = 3, scores = 'regression') #укажем scores, чтобы сохранились сами значения факторов
fan

## 
## Call:
## factanal(x = pizza2, factors = 3, scores = "regression")
## 
## Uniquenesses:
##   mois   prot    fat    ash sodium   carb    cal 
##  0.005  0.005  0.005  0.091  0.114  0.005  0.005 
## 
## Loadings:
##        Factor1 Factor2 Factor3
## mois            0.266  -0.961 
## prot    0.199   0.971  -0.116 
## fat     0.904   0.356   0.232 
## ash     0.629   0.706  -0.118 
## sodium  0.889   0.278   0.135 
## carb   -0.539  -0.714   0.443 
## cal     0.632           0.772 
## 
##                Factor1 Factor2 Factor3
## SS loadings      2.736   2.228   1.815
## Proportion Var   0.391   0.318   0.259
## Cumulative Var   0.391   0.709   0.968
## 
## Test of the hypothesis that 3 factors are sufficient.
## The chi square statistic is 2281.07 on 3 degrees of freedom.
## The p-value is 0

О, норм, работает.

Итак, что видно в аутпуте: уникальности, матрица факторных нагрузок, доля объясненной дисперсии и даже какой-то статистический тест с p-value.

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

Смотрим на таблицу факторных нагрузок. Так получается, что она еще и автоматически корреляционная матрица переменных с факторами. Вот почему:

\[ \mathrm{corr}(\boldsymbol{X}_i, \boldsymbol{F}_j) = \mathrm{corr}(a_{i1} \boldsymbol{F}_1 + a_{i2} \boldsymbol{F}_2 + \dots + a_{ip} \boldsymbol{F}_p + \boldsymbol{U}_i, \boldsymbol{F}_j) = \\ = \mathrm{corr}(a_{i1} \boldsymbol{F}_1, \boldsymbol{F}_j) + \mathrm{corr}(a_{i2} \boldsymbol{F}_2, \boldsymbol{F}_j) + \dots + \mathrm{corr}(a_{ij}\boldsymbol{F}_j, \boldsymbol{F}_j) + \dots + \mathrm{corr}(a_{ip} \boldsymbol{F}_p, \boldsymbol{F}_j) + \mathrm{corr}(\boldsymbol{U}_i, \boldsymbol{F}_j) = \\ = a_{ij} \, \mathrm{corr}(\boldsymbol{F}_j, \boldsymbol{F}_j) = a_{ij} \, \mathrm{var}(\boldsymbol{F}_j) = a_{ij} \]

Таким образом, мы интерпретируем факторы так же, как мы интерпретировали главные компоненты — на основании корреляций с исходными переменными. Как бы вы проинтерпретировали получившиеся три фактора?

Можно подобрать модель с двумя факторами:

factanal(pizza2, factors = 2)

## 
## Call:
## factanal(x = pizza2, factors = 2)
## 
## Uniquenesses:
##   mois   prot    fat    ash sodium   carb    cal 
##  0.005  0.520  0.005  0.191  0.126  0.071  0.005 
## 
## Loadings:
##        Factor1 Factor2
## mois    0.406  -0.911 
## prot    0.684  -0.111 
## fat     0.827   0.558 
## ash     0.895         
## sodium  0.807   0.472 
## carb   -0.934   0.238 
## cal     0.274   0.960 
## 
##                Factor1 Factor2
## SS loadings      3.717   2.364
## Proportion Var   0.531   0.338
## Cumulative Var   0.531   0.869
## 
## Test of the hypothesis that 2 factors are sufficient.
## The chi square statistic is 3485.54 on 8 degrees of freedom.
## The p-value is 0

Как теперь их можно интерпретировать? Отличается ли результат от PCA?

Из результатов факторного анализа можно отдельно доставать всякие штуки:

fan$loadings # факторные нагрузки aka матрица корреляций переменных с факторами

## 
## Loadings:
##        Factor1 Factor2 Factor3
## mois            0.266  -0.961 
## prot    0.199   0.971  -0.116 
## fat     0.904   0.356   0.232 
## ash     0.629   0.706  -0.118 
## sodium  0.889   0.278   0.135 
## carb   -0.539  -0.714   0.443 
## cal     0.632           0.772 
## 
##                Factor1 Factor2 Factor3
## SS loadings      2.736   2.228   1.815
## Proportion Var   0.391   0.318   0.259
## Cumulative Var   0.391   0.709   0.968

fan$correlation # корреляционная матрица исходных переменных

##              mois       prot        fat        ash     sodium        carb
## mois    1.0000000  0.3602477 -0.1713182  0.2655555 -0.1022789 -0.59180165
## prot    0.3602477  1.0000000  0.4980017  0.8238437  0.4291295 -0.85354226
## fat    -0.1713182  0.4980017  1.0000000  0.7916340  0.9333252 -0.64023817
## ash     0.2655555  0.8238437  0.7916340  1.0000000  0.8081221 -0.89898837
## sodium -0.1022789  0.4291295  0.9333252  0.8081221  1.0000000 -0.62017634
## carb   -0.5918017 -0.8535423 -0.6402382 -0.8989884 -0.6201763  1.00000000
## cal    -0.7644405  0.0702581  0.7645671  0.3264685  0.6719575 -0.02348458
##                cal
## mois   -0.76444054
## prot    0.07025810
## fat     0.76456710
## ash     0.32646845
## sodium  0.67195750
## carb   -0.02348458
## cal     1.00000000

fan$uniquenesses # уникальности

##       mois       prot        fat        ash     sodium       carb 
## 0.00500000 0.00500000 0.00500000 0.09123351 0.11363405 0.00500000 
##        cal 
## 0.00500000

fan$factors # количество факторов

## [1] 3

fan$dof # число степеней сводобы

## [1] 3

fan$rotmat # rotation matrix (кстати, мы по умолчанию применили varimax)

##            [,1]      [,2]       [,3]
## [1,]  0.7536936 0.1565949  0.6382978
## [2,]  0.1916326 0.8766368 -0.4413444
## [3,] -0.6286676 0.4549571  0.6307068

fan$n.obs # число наблюдений

## [1] 300

fan$scores # значения факторов

##              Factor1     Factor2      Factor3
##   [1,]  2.2403606726  1.01365654  1.586355336
##   [2,]  2.1517548542  1.00091815  1.496408251
##   [3,]  2.4850667309  0.72364156  1.437571756
##   [4,]  2.1961101438  0.77695944  1.231803288
##   [5,]  1.9063377764  1.02361718  1.317077547
##   [6,]  2.1080650467  0.80291517  1.169719507
##   [7,]  1.9106978395  0.96219229  1.226313169
##   [8,]  1.8419336983  1.08904918  1.501611777
##   [9,]  2.3597231318  0.82839938  1.490011356
##  [10,]  2.2997558996  0.96763399  1.580551579
##  [11,]  2.3171874472  0.96533982  1.625380205
##  [12,]  2.2394984222  0.99187466  1.670151354
##  [13,]  2.4010735791  0.73763937  1.413621437
##  [14,]  2.1691001728  0.69065729  1.192213382
##  [15,]  2.0138510469  0.91483867  1.234337039
##  [16,]  2.5802939761  0.32831412  1.032800530
##  [17,]  1.8510266753  0.81639109  0.990412118
##  [18,]  2.2745608639  0.28807009  0.632433234
##  [19,]  1.8753957387  0.87708898  1.099961226
##  [20,]  2.2117838012  0.59714566  1.099454189
##  [21,]  2.6617777924  0.62981107  1.556565161
##  [22,]  2.4447180621  0.73872153  1.508952001
##  [23,]  2.0620657018  0.87473489  1.221283371
##  [24,]  2.3612603970  0.57294296  1.001523050
##  [25,]  2.3244039216  0.38111229  0.802152569
##  [26,]  2.3662367303  0.68194554  1.214699229
##  [27,]  2.4256829077  0.72226741  1.389241102
##  [28,]  2.4316853463  0.49596114  1.155773681
##  [29,]  2.1243633904  0.96229702  1.377906129
##  [30,]  1.5245191852 -0.38591847 -1.115422504
##  [31,]  1.1144121194 -0.23128326 -1.381128600
##  [32,]  1.6946707363 -0.49788000 -1.033852276
##  [33,]  1.4371994153 -0.36068816 -1.171080576
##  [34,]  1.4071001078 -0.45100948 -1.255734239
##  [35,]  1.5788200815 -0.46090330 -1.181013939
##  [36,]  1.7296881919 -0.56912880 -1.298799704
##  [37,]  1.7466188162 -0.64082886 -1.294966697
##  [38,]  1.1836832066 -0.44264630 -1.587410267
##  [39,]  1.0101788214 -0.20634894 -1.317818517
##  [40,]  1.5337258621 -0.35655449 -1.179478042
##  [41,]  1.5888661077 -0.41395683 -1.124934026
##  [42,]  1.3317229696 -0.31112092 -1.262173758
##  [43,]  1.5426171736 -0.46663798 -1.163324108
##  [44,]  0.8407773262 -0.22250609 -1.622066201
##  [45,]  0.6487965309 -0.08028165 -1.593806926
##  [46,]  1.0805401304 -0.29427053 -1.547582136
##  [47,]  1.5529898252 -0.37247076 -1.083277985
##  [48,]  1.5709499707 -0.49863089 -1.277396981
##  [49,]  0.9984879853 -0.23423108 -1.376616855
##  [50,]  1.5152058510 -0.43750266 -1.169501388
##  [51,]  1.6089133152 -0.48940577 -1.144758288
##  [52,]  1.6130042718 -0.50646770 -1.124521755
##  [53,]  1.3698388494 -0.34341836 -1.247357567
##  [54,]  1.3595486121 -0.33410616 -1.261472695
##  [55,]  1.0559521732 -0.28792381 -1.545025210
##  [56,]  1.3435712796 -0.38944131 -1.400085217
##  [57,]  1.3979053333 -0.34587471 -1.338457029
##  [58,]  1.4837380311 -0.39847340 -1.329807264
##  [59,]  1.5837908950 -0.52816023 -1.371221916
##  [60,]  1.5254965452 -0.50539511 -1.356971003
##  [61,] -0.6564304681  2.11766623 -0.166976218
##  [62,] -1.3955566660  2.67209895 -0.024385555
##  [63,] -0.6439029608  1.97154837 -0.256448225
##  [64,] -0.8744255312  2.08895038 -0.239432061
##  [65,] -0.2090670311  1.66833341 -0.343001080
##  [66,] -0.7699740622  2.03786028 -0.293637694
##  [67,] -1.5347278184  2.65119362 -0.190832803
##  [68,] -1.1917307459  2.55270925 -0.003793950
##  [69,] -1.3468090115  2.45350710 -0.181656518
##  [70,] -1.2062193773  2.22601950 -0.487925867
##  [71,] -1.4109000180  2.27302994 -0.535696204
##  [72,] -0.7438030121  1.98352299 -0.290824607
##  [73,] -0.6436995362  2.12107486 -0.131094138
##  [74,] -0.7390440384  2.10206079 -0.109591606
##  [75,] -0.6843265470  2.04250364 -0.117166835
##  [76,] -0.5667719396  1.83401273 -0.301262709
##  [77,] -0.6512577379  1.77796905 -0.441735958
##  [78,] -1.1236981168  2.17386524 -0.593206084
##  [79,] -0.7317686173  2.16958936 -0.006643861
##  [80,] -0.9075375101  2.29165938 -0.132728822
##  [81,] -1.1963530935  2.23114877 -0.561171329
##  [82,] -1.1605763707  2.32554134 -0.167078126
##  [83,] -1.0233565303  1.97503662 -0.546891922
##  [84,] -1.2919992535  2.46241444 -0.128349471
##  [85,] -0.1391120436  1.67880770 -0.190705740
##  [86,] -1.2489448170  2.25579590 -0.386383866
##  [87,] -0.6611091064  2.00245912 -0.307434785
##  [88,]  0.1885979781  1.21958461 -0.326288126
##  [89,]  0.2601249140  1.24838692 -0.234436763
##  [90,] -0.3538261531  1.42592007 -0.337197495
##  [91,] -0.6017589423  1.54298063 -0.299575284
##  [92,] -0.3459991818  1.38668573 -0.275651040
##  [93,] -0.0146485657  1.31432242 -0.322070333
##  [94,] -0.3773789287  1.38084025 -0.336722579
##  [95,] -0.3654960133  1.40428407 -0.388891708
##  [96,] -1.3255604922  2.56714064 -0.094706030
##  [97,] -0.3613116928  1.46871134 -0.275505996
##  [98,] -0.2476722068  1.31202060 -0.338505886
##  [99,] -0.2204720598  1.24212636 -0.383045628
## [100,] -0.5767082151  1.50063565 -0.338952718
## [101,] -0.2882175965  1.51884367 -0.417630932
## [102,] -0.3010984843  1.48050831 -0.394382685
## [103,] -0.3476628457  1.53855072 -0.449358034
## [104,] -0.3809592875  1.44427420 -0.313896243
## [105,] -0.4242264947  1.40291836 -0.353980844
## [106,] -0.4148148445  1.42966395 -0.385380757
## [107,] -0.2849448624  1.28068192 -0.426206521
## [108,] -0.4076482165  1.48865177 -0.299168407
## [109,] -0.4703742257  1.49515330 -0.331947871
## [110,] -0.5837414082  1.54017359 -0.378226432
## [111,]  0.2228909974  1.18798010 -0.289467744
## [112,] -0.6630063572  1.53609572 -0.431809952
## [113,] -0.2491017302  1.55139606 -0.495807157
## [114,] -0.4699016029  1.47621685 -0.351585932
## [115,]  0.2094440488  1.24571710 -0.274406306
## [116,]  0.3007729562  1.21284206 -0.295522536
## [117,]  0.3331837834  1.09597392 -0.404659545
## [118,]  0.2616061627  1.09241178 -0.381803577
## [119,] -1.1236981168  2.17386524 -0.593206084
## [120,] -0.3000530023 -0.83158171  0.488023463
## [121,] -0.1318964137 -0.84396801  0.237619511
## [122,] -0.3913683952 -0.71954937  0.459104420
## [123,] -1.7824889326 -0.37073958  0.120606209
## [124,] -0.5062646322 -0.73396282  0.521536476
## [125,] -1.7508036890 -0.43167666  0.069029838
## [126,] -0.6599641614 -0.69749429  0.487625793
## [127,] -0.1471450550 -0.82058298 -0.033604901
## [128,] -0.3672733856 -0.78693893  0.332512908
## [129,] -0.3268926424 -0.80095526  0.275289584
## [130,] -0.6628668564 -0.65174582  0.184293672
## [131,] -0.6314139563 -0.68847351  0.446553197
## [132,]  0.0620842391 -0.87847252  0.476208521
## [133,]  0.0866662887 -0.88267477  0.614038597
## [134,]  0.2646732156 -1.00485599  0.115112698
## [135,] -0.5758860162 -0.73169183  0.188392424
## [136,] -0.7291339553 -0.74901749  0.294293682
## [137,]  0.2646732156 -1.00485599  0.115112698
## [138,] -0.5689581078 -0.68164815  0.556653854
## [139,] -0.4412783353 -0.76979556  0.268819563
## [140,] -0.1092327637 -0.81326762  0.511861391
## [141,] -0.0714279440 -0.88920515  0.353315373
## [142,]  0.2705697644 -1.04010003  0.347424840
## [143,]  0.1322437936 -0.96480380  0.533364254
## [144,]  0.0690241268 -0.93870178  0.418985182
## [145,] -1.2472538982 -0.66558367  0.228233689
## [146,] -0.6628668564 -0.65174582  0.184293672
## [147,] -0.6131405536 -0.74431085  0.518370420
## [148,] -0.2633341239 -0.71595766  1.217767729
## [149,] -0.6452732493 -0.60657008  1.035891625
## [150,] -0.7863210938 -0.54181156  1.082516215
## [151,]  0.0399793183 -0.75014649  0.989664128
## [152,] -0.1958271874 -0.72955635  1.172809228
## [153,] -0.6889072014 -0.58321943  1.017012897
## [154,] -0.5534773902 -0.65909045  1.172386875
## [155,] -0.8796835042 -0.51342908  0.968200779
## [156,] -0.7625457389 -0.54292645  1.082984946
## [157,] -0.5136323523 -0.67488901  0.992410797
## [158,] -0.1260441128 -0.69718469  1.016658593
## [159,] -0.6630200341 -0.62791921  1.091423664
## [160,] -0.1433731874 -0.67466674  1.133934202
## [161,] -0.6633503098 -0.60652431  1.055744228
## [162,] -0.4456240518 -0.65538122  1.302235507
## [163,] -0.3299463103 -0.67951831  1.202957460
## [164,] -0.0462946668 -0.71034224  1.106137789
## [165,] -0.7922443006 -0.54985709  1.111121538
## [166,] -0.2765030338 -0.66555894  1.392900093
## [167,] -0.7643990265 -0.58073080  0.994930751
## [168,] -0.9308335754 -0.53026450  0.992796345
## [169,] -1.0378602881 -0.50927628  1.060276180
## [170,] -0.5252360474 -0.61711983  1.118287551
## [171,] -0.5604205081 -0.67371934  1.046086695
## [172,] -0.3016177534 -0.74021680  1.105263592
## [173,] -0.4620008841 -0.67775527  1.132068625
## [174,] -0.4152371410 -0.65820714  1.153426690
## [175,] -0.2392387783 -0.73462895  1.178923743
## [176,] -0.1347829908 -0.75465222  1.170848548
## [177,] -0.7065606368 -0.66007230  0.907077888
## [178,] -0.7052084389 -0.52907926  1.270454126
## [179,] -0.7669398263 -0.52418763  1.254392851
## [180,] -0.8287004468 -0.56989612  0.965871034
## [181,] -0.7682631495 -0.50400407  1.275308271
## [182,] -0.0681737908 -0.67197746  1.153811750
## [183,] -0.5204759038 -0.60207280  1.235604840
## [184,] -0.6092084547 -0.54575639  1.207813338
## [185,] -0.6972704977 -0.51567524  1.314716201
## [186,] -0.5939143569 -0.54508464  1.278438987
## [187,] -0.6157929001 -0.53015748  1.272934448
## [188,] -0.9075546254 -0.45209687  1.298956859
## [189,] -0.9034021373 -0.42248813  1.284882425
## [190,] -0.7331520605 -0.49758105  1.228558551
## [191,] -0.9381055131 -0.60315974  0.720210141
## [192,] -0.1542219603 -0.65658458  1.424124840
## [193,] -0.7970942874 -0.44354121  1.388287278
## [194,] -0.8344037197 -0.48926991  1.343928463
## [195,] -0.7272484866 -0.54664030  1.312450364
## [196,] -0.3634201019 -0.63934454  0.970983776
## [197,] -0.8669825809 -0.46824072  1.186391474
## [198,] -0.4458325684 -0.58687218  1.202331935
## [199,] -0.8403950759 -0.47135000  1.380179771
## [200,] -0.9311263837 -0.49171184  1.261024081
## [201,] -0.4837883396 -0.58944642  1.460393423
## [202,] -0.4120952211 -0.62578548  1.429259752
## [203,] -0.7603225496 -0.50858031  1.320757323
## [204,] -0.6511242527 -0.45611815  1.633198829
## [205,] -1.1083210972 -0.41517655  1.198410938
## [206,] -0.5915123496 -0.59547625  1.196299077
## [207,] -0.8966051013 -0.57599698  0.396145385
## [208,]  0.0466019595 -0.81639803  0.433541167
## [209,] -0.3893761776 -0.78932624  0.400991460
## [210,] -1.7164518588 -0.27511957  0.544040743
## [211,] -0.3265123805 -0.73624683  0.403224130
## [212,] -0.2932551471 -0.75588575  0.392583499
## [213,] -1.6418119782 -0.47483164 -0.065190853
## [214,] -0.2319282837 -0.86766588  0.350593553
## [215,] -0.4093341372 -0.80398693  0.656136555
## [216,] -0.5611281783 -0.66226519  0.627610390
## [217,] -1.8047530753 -0.38026537  0.170378061
## [218,] -0.0455772260 -0.86881185  0.591489371
## [219,] -0.4972289904 -0.61334498  0.354144355
## [220,] -0.5771739885 -0.71612835  0.389513271
## [221,] -0.7870893792 -0.60044211  0.337399891
## [222,] -0.6904069626 -0.65067737  0.444019329
## [223,] -0.4688860134 -0.76426230  0.259689750
## [224,] -0.3526204332 -0.81986949  0.382351195
## [225,] -0.4063563585 -0.74899600  0.583501109
## [226,] -0.5733465959 -0.72541362  0.380232698
## [227,] -0.5483667556 -0.72460286  0.391720562
## [228,]  0.1224767974 -0.93423499  0.323200478
## [229,] -0.2640182510 -0.73059744  0.358039524
## [230,] -0.6152917257 -0.70061070  0.426118720
## [231,] -0.2207262896 -0.86301422  0.308637907
## [232,] -0.3533345249 -0.79725501  0.336647802
## [233,]  0.0316932063 -0.93931308  0.224401367
## [234,] -0.4561786417 -0.80227378  0.338105783
## [235,] -0.7760109052 -0.65134221  0.354519098
## [236,] -0.0584421152 -0.92914021  0.452832672
## [237,] -1.1641433622 -0.60578215  0.316331489
## [238,] -0.6548454950 -0.68458555  0.380456991
## [239,] -0.5483667556 -0.72460286  0.391720562
## [240,] -0.2261682399 -0.62774049 -1.657090791
## [241,] -0.1053217215 -0.73691400 -1.666798111
## [242,] -0.2273248664 -0.64091000 -1.564983928
## [243,] -0.1193538843 -0.60003667 -1.602905105
## [244,] -0.0891131708 -0.66492303 -1.472197597
## [245,] -0.4969962743 -0.45410882 -1.632509013
## [246,] -0.2483464369 -0.85944509 -1.887226594
## [247,] -0.0970237888 -0.74360357 -1.858042502
## [248,] -0.2641877399 -0.73832715 -1.663799441
## [249,] -0.0049763509 -0.55948273 -1.538010395
## [250,] -0.1125711670 -0.65963883 -1.671613523
## [251,] -0.2310873737 -0.47450775 -1.444859445
## [252,] -0.3081643611 -0.52574339 -1.679093153
## [253,] -0.1593471933 -0.57962511 -1.662364206
## [254,] -0.2789489329 -0.61110059 -1.736717879
## [255,] -0.0087151472 -0.79989235 -1.804803439
## [256,] -0.2258723792 -0.69881743 -1.658660575
## [257,] -0.2512359869 -0.65453256 -1.666814956
## [258,] -0.3690831334 -0.54182691 -1.570913580
## [259,] -0.3274008558 -0.51338581 -1.413770692
## [260,] -0.3064494648 -0.43113070 -1.564645325
## [261,] -0.1701732230 -0.57316168 -1.616412981
## [262,] -0.1945782780 -0.62726961 -1.637032756
## [263,] -0.2005364068 -0.57558560 -1.569821295
## [264,] -0.3343848744 -0.75388906 -1.963691419
## [265,] -0.1700467032 -0.60998592 -1.659556801
## [266,] -0.1865690767 -0.66843038 -1.612812090
## [267,] -0.3176273572 -0.53318695 -1.507685396
## [268,] -0.1586424822 -0.73766474 -1.755852980
## [269,] -0.0912035483 -0.54523205 -0.733212692
## [270,] -0.1896130982 -0.38994722 -0.429646993
## [271,] -0.1864072779 -0.52926831 -0.873944328
## [272,] -0.1264472341 -0.61768326 -0.817252368
## [273,] -0.0914707391 -0.62231127 -0.732443888
## [274,] -0.2400717376 -0.56506499 -0.909002089
## [275,]  0.0023752220 -0.53682294 -0.687442611
## [276,] -0.0450163720 -0.49768748 -0.749746599
## [277,]  0.0139827070 -0.50856352 -0.712657273
## [278,]  0.0248556056 -0.69824580 -1.029945432
## [279,] -0.4482304001 -0.03404786 -0.412792955
## [280,] -0.4227670694 -0.39271346 -0.674079591
## [281,] -0.0145035502 -0.56212568 -0.869975519
## [282,] -0.1981420653 -0.37251769 -0.682075736
## [283,] -0.1319332508 -0.46869511 -0.603371462
## [284,] -0.2032093971 -0.51540468 -0.812895514
## [285,]  0.2166337669 -0.49911115 -0.426480989
## [286,]  0.0027125211 -0.60372788 -0.920509837
## [287,] -0.0009888876 -0.53971526 -0.616514538
## [288,] -0.1600400524 -0.48270054 -0.559359510
## [289,] -0.0815784835 -0.50492781 -0.902831050
## [290,] -0.0587008082 -0.56497971 -0.675739390
## [291,] -0.0643984553 -0.55783784 -0.725236942
## [292,] -0.1443051554 -0.55601196 -0.859305404
## [293,] -0.0477442039 -0.69937331 -0.659514402
## [294,]  0.1458822477 -0.65520493 -0.778702613
## [295,]  0.0751389883 -0.45088887 -0.484026429
## [296,] -0.1249121426 -0.41300520 -0.564966576
## [297,] -0.0184274013 -0.28697552 -0.324125970
## [298,] -0.2266289174 -0.38575999 -0.506405446
## [299,] -0.1864072779 -0.52926831 -0.873944328
## [300,] -0.1264472341 -0.61768326 -0.817252368

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

library(plotly)

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

plot_ly(x = fan$scores[,1], y = fan$scores[,2], z = fan$scores[,3], color = pizza$brand,
                type = 'scatter3d', mode = 'markers') %>% 
  layout(scene = list(
    xaxis = list(title = 'Factor 1'),
    yaxis = list(title = 'Factor 2'),
    zaxis = list(title = 'Factor 3')
  ))

## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors

## Warning in RColorBrewer::brewer.pal(N, "Set2"): n too large, allowed maximum for palette Set2 is 8
## Returning the palette you asked for with that many colors