PCA: Как работает princomp() и могу ли я использовать его для выбора переменных для ARIMA?

Question

Я пытаюсь использовать PCA, чтобы выбрать хороших предикторов для использования в xreg аргумент arima модель, чтобы попытаться предсказать tVar переменная ниже. Я просто использую приведенный ниже набор данных с несколькими переменными, чтобы сделать пример простым.

Я пытаюсь понять, как формула аргумента в princomp работает. Для pc объект ниже, это говорит "использовать xVar1 а также xVar2 объяснить разницу в na.omit(dfData[,c("tVar","xVar1","xVar2")])"?

В конечном итоге я хотел бы создать новую переменную, которая объясняет большую часть различий в tVar, Это то, что я могу сделать с помощью PCA? Если да, может кто-нибудь объяснить, как или указать мне на пример?

Код:

pc <- princomp(~xVar1+xVar2,
               data = na.omit(dfData[,c("tVar","xVar1","xVar2")]), 
               cor=TRUE)

Данные:

dput(na.omit(dfData[1:100,c("tVar","xVar1","xVar2")]))
structure(list(tVar = c(11, 14, 17, 5, 5, 5.5, 8, 5.5, 
          6.5, 8.5, 4, 5, 9, 10, 11, 7, 6, 7, 7, 5, 6, 9, 9, 6.5, 9, 3.5, 
          2, 15, 2.5, 17, 5, 5.5, 7, 6, 3.5, 6, 9.5, 5, 7, 4, 5, 4, 9.5, 
          3.5, 5, 4, 4, 9, 4.5, 6, 10, 9.5, 15, 9, 5.5, 7.5, 12, 17.5, 
          19, 7, 14, 17, 3.5, 6, 15, 11, 10.5, 11, 13, 9.5, 9, 7, 4, 6, 
          15, 5, 18, 5, 6, 19, 19, 6, 7, 7.5, 7.5, 7, 6.5, 9, 10, 5.5, 
          5, 7.5, 5, 4, 10, 7, 5, 12), xVar1 = c(0L, 0L, 0L, 0L, 
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 
          1L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L),
          xVar2  = c(0L, 
          1L, 0L, 1L, 0L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 1L, 0L, 1L, 1L, 0L, 
          2L, 3L, 0L, 0L, 1L, 0L, 0L, 1L, 1L, 0L, 1L, 0L, 0L, 0L, 1L, 0L, 
          0L, 1L, 0L, 1L, 1L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 0L, 1L, 0L, 0L, 
          0L, 0L, 0L, 0L, 1L, 0L, 1L, 1L, 0L, 0L, 0L, 3L, 1L, 0L, 1L, 2L,
          0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 0L, 1L, 1L, 1L, 0L, 
          1L, 1L, 0L, 0L, 1L, 0L, 1L, 1L, 0L, 0L, 1L, 0L, 0L, 0L, 0L, 1L, 
          0L)), .Names = c("tVar", "xVar1", "xVar2"
          ), row.names = c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 9L, 10L, 11L, 12L, 
          13L, 14L, 15L, 16L, 17L, 18L, 19L, 20L, 21L, 22L, 23L, 24L,25L, 
          26L, 27L, 28L, 29L, 30L, 31L, 32L, 33L, 34L, 35L, 36L, 37L,38L, 
          39L, 40L, 41L, 42L, 43L, 44L, 45L, 46L, 47L, 48L, 49L, 50L,51L, 
          52L, 54L, 55L, 56L, 57L, 58L, 59L, 60L, 61L, 62L, 63L, 64L, 65L,
          66L, 67L, 68L, 69L, 70L, 71L, 72L, 73L, 74L, 75L, 76L, 77L, 78L, 
          79L, 80L, 81L, 82L, 83L, 84L, 85L, 86L, 87L, 88L, 89L, 90L, 91L, 
          92L, 93L, 94L, 95L, 96L, 97L, 98L, 99L, 100L),
          class  = "data.frame", na.action = structure(c(8L,53L),
          .Names = c("8", "53"), class = "omit"))

Answer 1

(Это очень хороший пост. Интересно иметь еще один пост сегодня о PCA. Хотя этот вопрос является более основным, относительно разницы между princomp а также prcomp , но математические детали с кодом R, который я делаю в ответе, могут быть полезны для любого изучающего PCA.)

PCA используется для уменьшения размера (приближение низкого ранга), когда:

1. у тебя много (скажем p) коррелированные переменные x1, x2, ..., xp;
2. Вы хотите сократить их до небольшого числа (скажем, k < p) новых, линейно независимых переменных z1, z2, ..., zk;
3. ты хочешь использовать z1, z2, ..., zk скорее, чем x1, x2, ..., xp предсказать переменную ответа y,

---

Фундаментальная картина и немного математики

Предположим, у вас есть переменная ответа y, полная линейная регрессия без отбрасывания каких-либо переменных должна принимать формулу:

y ~ x1 + x2 + ... + xp

Тем не менее, мы можем сделать разумную приблизительную модель, после PCA. Позволять X быть матрицей модели выше, т. е. матрицей путем объединения всех наблюдений x1, x2, ... , xp по столбцу, то

S <- cor(X)  ## get correlation matrix S
E <- eigen(S)  ## compute eigen decomposition of S
root_eigen_value <- sqrt(E$values)  ## square root of eigen values
eigen_vector_mat <- E$vectors  ## matrix of eigen vectors
X1 <- scale(X) %*% eigen_vector_mat  ## transform original matrix

Сейчас, root_eigen_value (длина- p vector) монотонно убывает, т. е. вклад в полную ковариацию уменьшается, поэтому мы можем выбрать только первое k ценности. Соответственно, мы можем выбрать первый k столбцы преобразованной матрицы X1, Давайте сделаем:

Z <- X1[, 1:k]

Теперь мы успешно сократили p переменные к k переменные, и каждый столбец Z это новая переменная z1, z2, ..., zk, Имейте в виду, что эти переменные не являются подмножеством исходных переменных; они совершенно новые, без имен. Но так как мы заинтересованы только в прогнозировании y, не имеет значения, какое имя мы даем z1, z2, ..., zk, Тогда мы можем подобрать приближенную линейную модель:

y ~ z1 + z2 + ... + zk

---

использование princomp()

На самом деле все проще, потому что princomp() делает все вычисления для нас. По телефону:

pc <- princomp(~ x1 + x2 + ... + xp, data, cor = TRUE)

мы можем получить все, что хотим. Среди нескольких возвращаемых значений в pc:

1. pc$sdev дает root_eigen_value, Если вы делаете plot(pc) Вы можете увидеть барплот, показывающий это. Если ваши входные данные сильно коррелированы, то на этом рисунке ожидается, что вы увидите почти экспоненциальный спад, причем только несколько переменных доминируют над ковариацией. (К сожалению, ваши данные игрушки не будут работать. xVar1 а также xVar2 являются двоичными, и они уже линейно независимы, следовательно, после PCA, вы увидите, что они оба дают равный вклад.)
2. pc$loadings дает eigen_vector_mat;
3. pc$scores дает X1,

---

использование arima()

Процесс выбора переменных прост. Если вы решили взять первый k переменные из общего числа p переменные, проверяя plot(pc) затем вы извлекаете первый k колонны pc$scores матрица. Каждый столбец формы z1, z2, ..., zk и передать их arima() через аргумент reg,

---

Вернуться к вашему вопросу о формуле

Для объекта pc, приведенного ниже, он говорит "используйте xVar1 и xVar2, чтобы объяснить разницу в na.omit(dfData[,c("tVar","xVar1","xVar2")])"

После моего объяснения вы должны знать, что ответ "Нет". Не смешивайте переменную ответа tVar используется на шаге регрессии с переменными предиктора xVar1, xVars... используется на этапе PCA.

princomp() позволяет передавать аргументы тремя способами:

1. по формуле и данным;
2. по матрице моделей;
3. по ковариационной матрице.

Вы выбрали первый путь. Формула используется, чтобы сказать princomp() извлечь данные из data и позже он будет вычислять матрицу модели, ковариационную матрицу, корреляционную матрицу, собственное разложение, пока мы, наконец, не получим результат PCA.

---

Отслеживание ваших комментариев

Так что, если я правильно понимаю, PCA в первую очередь для уменьшения количества переменных, и я не должен включать переменную ответа tVar в формуле или данных. Но мне было интересно, почему princomp(~xVar1+xVar2, data = na.omit(dfData[,c("tVar","xVar1","xVar2")]), cor=TRUE) а также princomp(na.omit(dfData[,c("xVar1","xVar2")]), cor=TRUE) в основном эквивалентны?

Формула рассказывает, как извлечь матрицу из фрейма данных. Так как вы используете ту же формулу ~ xVar1 + xVar2 включаете ли вы tVars в кадре данных для передачи в princomp не имеет значения, так как этот столбец не будет затронут princomp,

Не включать tVars в вашей формуле для PCA. Как я уже сказал, регрессия и PCA - это разные проблемы, и их не следует путать друг с другом.

Чтобы быть ясным, стратегия с PCA не заключается в создании новой переменной, которая является комбинацией xVar1 а также xVar2 и объясняет большинство различий в tVar, а скорее, чтобы создать новую переменную, которая является комбинацией xVar1 а также xVar2 и объясняет большинство различий dfData[,c("xVar1","xVar2")]?

Да. Регрессия (или arima() в ваших настройках) используется для установки связи между вашим ответом tVars и предикторные переменные x1, x2, ..., xp или же z1, z2, ..., zk, Модель регрессии / аримы объяснит среднее значение и дисперсию ответа в терминах предикторов.

Спс это другая проблема. Он только выбирает низкое ранговое (меньшее количество параметров) представление ваших исходных переменных предиктора xVar1, xVar2, ..., так что вы можете использовать меньше переменных в более поздней регрессии / моделировании ARIMA.

Тем не менее, вам, возможно, придется подумать, стоит ли вам делать PCA для решения вашей проблемы.

1. У вас много переменных, скажем, 10+? В статистическом моделировании обычно достигают сотен тысяч параметров. Вычисления могут стать очень медленными, если мы используем их все. PCA полезен в этом случае, чтобы уменьшить вычислительную сложность, в то же время давая разумное представление исходной ковариации.
2. Ваши переменные сильно коррелированы? Если они легко линейно независимы друг от друга, PCA не может ничего отбросить. Например, данные игрушки xVar1 а также xVar2 Вы дали просто линейно независимы, поэтому уменьшение размеров невозможно. Вы можете просмотреть соотношение в ваших данных, pairs(mydata), Лучшей визуализацией может быть использование corrplot R пакет. Посмотрите этот ответ для примеров того, как использовать его для построения ковариационной матрицы.