Тестирование на линейную зависимость среди столбцов матрицы

26

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

Моя гипотеза состоит в том, что по крайней мере одна ценная бумага линейно зависит от других ценных бумаг. Есть ли в R функция, которая последовательно проверяет каждый столбец на наличие матрицы линейной зависимости?

Например, один из подходов заключается в построении корреляционной матрицы по одной безопасности за раз и вычислении детерминанта на каждом этапе. Когда определитель = 0, остановитесь, поскольку вы определили ценную бумагу, которая представляет собой линейную комбинацию других ценных бумаг.

Любые другие методы для определения линейной зависимости в такой матрице приветствуются.

r correlation pca linear-model svd Рам Ахлувалия
источник

Ваша матрица является положительно-полуопределенной, но не положительно-определенной, поскольку она является единственной.

ttnphns

Каковы размеры (нет. Переменных; нет. Образцов)?

Карл

Количество столбцов = 480. Количество строк для каждого временного ряда = 502. В целом вы обнаружите, что чем больше временной ряд, тем ковариационная матрица выборки имеет тенденцию быть положительно определенной. Однако во многих случаях вы хотите использовать существенно меньшее значение T (или экспоненциальный вес), чтобы отразить последние рыночные условия.

Рам Ахлувалия

3

Вопрос некорректен. Если ваша матрица данных равна 480 на 502, то сказать, что матрица имеет ранг

(пространство столбцов матрицы имеет размерность

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

q < 480

$q < 480$

q < 480

$q < 480$

NRH

Ковариационная матрица симметрична. Он генерируется транспонированием (A) * A. Матрица A имеет размеры 480x502. Однако ковариационная матрица имеет размер 480x480

Рам Ахлувалия,

6

Вы, кажется, задаете действительно провокационный вопрос: как определить, учитывая матрицу единственной корреляции (или ковариации, или суммы квадратов и перекрестных произведений), какой столбец линейно зависит от какого. Я предполагаю, что операция очистки может помочь. Вот мое исследование в SPSS (не R), чтобы проиллюстрировать.

Давайте сгенерируем некоторые данные:

        v1        v2        v3         v4          v5
    -1.64454    .35119   -.06384    -1.05188     .25192
    -1.78520   -.21598   1.20315      .40267    1.14790
     1.36357   -.96107   -.46651      .92889   -1.38072
     -.31455   -.74937   1.17505     1.27623   -1.04640
     -.31795    .85860    .10061      .00145     .39644
     -.97010    .19129   2.43890     -.83642    -.13250
     -.66439    .29267   1.20405      .90068   -1.78066
      .87025   -.89018   -.99386    -1.80001     .42768
    -1.96219   -.27535    .58754      .34556     .12587
    -1.03638   -.24645   -.11083      .07013    -.84446

Давайте создадим некоторую линейную зависимость между V2, V4 и V5:

compute V4 = .4*V2+1.2*V5.
execute.

Итак, мы изменили нашу колонку V4.

matrix.
get X. /*take the data*/
compute M = sscp(X). /*SSCP matrix, X'X; it is singular*/
print rank(M). /*with rank 5-1=4, because there's 1 group of interdependent columns*/
loop i= 1 to 5. /*Start iterative sweep operation on M from column 1 to column 5*/
-compute M = sweep(M,i).
-print M. /*That's printout we want to trace*/
end loop.
end matrix.

Распечатки М в 5 итераций:

M
     .06660028    -.12645565    -.54275426    -.19692972    -.12195621
     .12645565    3.20350385    -.08946808    2.84946215    1.30671718
     .54275426    -.08946808    7.38023317   -3.51467361   -2.89907198
     .19692972    2.84946215   -3.51467361   13.88671851   10.62244471
     .12195621    1.30671718   -2.89907198   10.62244471    8.41646486

M
     .07159201     .03947417    -.54628594    -.08444957    -.07037464
     .03947417     .31215820    -.02792819     .88948298     .40790248
     .54628594     .02792819    7.37773449   -3.43509328   -2.86257773
     .08444957    -.88948298   -3.43509328   11.35217042    9.46014202
     .07037464    -.40790248   -2.86257773    9.46014202    7.88345168

M
    .112041875    .041542117    .074045215   -.338801789   -.282334825
    .041542117    .312263922    .003785470    .876479537    .397066281
    .074045215    .003785470    .135542964   -.465602725   -.388002270
    .338801789   -.876479537    .465602725   9.752781632   8.127318027
    .282334825   -.397066281    .388002270   8.127318027   6.772765022

M
   .1238115070   .0110941027   .0902197842   .0347389906   .0000000000
   .0110941027   .3910328733  -.0380581058  -.0898696977  -.3333333333
   .0902197842  -.0380581058   .1577710733   .0477405054   .0000000000
   .0347389906  -.0898696977   .0477405054   .1025348498   .8333333333
   .0000000000   .3333333333   .0000000000  -.8333333333   .0000000000

M
   .1238115070   .0110941027   .0902197842   .0347389906   .0000000000
   .0110941027   .3910328733  -.0380581058  -.0898696977   .0000000000
   .0902197842  -.0380581058   .1577710733   .0477405054   .0000000000
   .0347389906  -.0898696977   .0477405054   .1025348498   .0000000000
   .0000000000   .0000000000   .0000000000   .0000000000   .0000000000

Обратите внимание, что в итоге столбец 5 заполнился нулями. Это означает (насколько я понимаю), что V5 линейно связан с некоторыми из предыдущих столбцов. Какие колонки? Посмотрите на итерацию, где последний столбец 5 не полон нулей - итерация 4. Мы видим там, что V5 связан с V2 и V4 с коэффициентами -.3333 и .8333: V5 = -.3333 * V2 + .8333 * V4, что соответствует к тому, что мы сделали с данными: V4 = .4 * V2 + 1.2 * V5.

Вот как мы узнали, какая колонка линейно связана с другой. Я не проверял, насколько полезен вышеуказанный подход в более общем случае со многими группами взаимозависимостей в данных. В приведенном выше примере это оказалось полезным, хотя.

ttnphns
источник

Разве это не уменьшенный рядный ряд? Если да, нет ли пакетов / функций, доступных в R?

Арун

@Arun, я не пользователь R, поэтому не могу знать.

ttnphns

25

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

В R:

rankifremoved <- sapply(1:ncol(your.matrix), function (x) qr(your.matrix[,-x])$rank)
which(rankifremoved == max(rankifremoved))

Джеймс
источник

1

Чрезвычайно полезный ответ при определении ошибочного столбца в матрице регрессии, где я получил ошибку system is exactly singular: U[5,5] = 0 , которая, как теперь я знаю, означает, что проблема возникла в столбце 5 (задним числом это кажется очевидным, поскольку это столбец нулей!)

Мэтт Веллер,

В комментарии Джеймса он разместил скрипт: rankifremoved <- sapply (1: ncol (your.matrix), function (x) qr (your.matrix [, - x]) $ rank), который (rankifremoved == max ( rankifremoved)) Я сделал тест на матрице, я хотел бы знать на выходе R. Столбцы вывода линейно зависимы? Поблагодарить!

@ EltonAraújo: на выходе будет вектор, дающий индексы линейно зависимых столбцов: так (2,4,5) для примера в ответе ttnphns. Но мне интересно, как вопросы числовой точности повлияют на этот метод.

Scortchi - Восстановить Монику

rankifremoved содержит все столбцы, которые линейно зависят между ними или между ними. В некоторых приложениях мы можем захотеть сохранить столбец или несколько столбцов и не отбрасывать все

MasterJedi

Разве это не должно возвращать пустой набор для your.matrix = matrix(1:4, 2)?

Хольгер Брандл

15

Вопрос заключается в том, чтобы «определить лежащие в основе [линейные] отношения» среди переменных.

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

$0$

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

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

pca <- function(x, threshold, ...) {
  fit <- princomp(x)
  #
  # Compute the relations among "small" components.
  #
  if(missing(threshold)) threshold <- max(fit$sdev) / ncol(x)
  i <- which(fit$sdev < threshold)
  relations <- fit$loadings[, i, drop=FALSE]
  relations <- round(t(t(relations) / apply(relations, 2, max)), digits=2)
  #
  # Plot the loadings, highlighting those for the small components.
  #
  matplot(x, pch=1, cex=.8, col="Gray", xlab="Observation", ylab="Value", ...)
  suppressWarnings(matplot(x %*% relations, pch=19, col="#e0404080", add=TRUE))

  return(t(relations))
}

$B,C,D,$ $E$ $A$

process <- function(z, beta, sd, ...) {
  x <- z %*% beta; colnames(x) <- "A"
  pca(cbind(x, z + rnorm(length(x), sd=sd)), ...)
}

$B, \ldots, E$ $A=B+C+D+E$ $A=B+(C+D)/2+E$ sweep

n.obs <- 80 # Number of cases
n.vars <- 4 # Number of independent variables
set.seed(17)
z <- matrix(rnorm(n.obs*(n.vars)), ncol=n.vars)
z.mean <- apply(z, 2, mean)
z <- sweep(z, 2, z.mean)
colnames(z) <- c("B","C","D","E") # Optional; modify to match `n.vars` in length

$B, \ldots, E$ $A$

Полученные результаты

Вывод, связанный с верхней левой панелью, был

       A  B  C  D  E
Comp.5 1 -1 -1 -1 -1

$0$ $0 \approx A -B-C-D-E$

Выход для верхней средней панели был

       A     B     C     D     E
Comp.5 1 -0.95 -1.03 -0.98 -1.02

$(A,B,C,D,E)$

       A     B     C     D     E
Comp.5 1 -1.33 -0.77 -0.74 -1.07

$A' = B' + C' + D' + E'$

$1, 1/2, 1/2, 1$

На практике это часто не тот случай, когда одна переменная выделяется как очевидная комбинация других: все коэффициенты могут иметь сопоставимые размеры и различные знаки. Более того, когда существует более одного измерения отношений, не существует уникального способа их указать: требуется дополнительный анализ (например, сокращение строк), чтобы определить полезную основу для этих отношений. Вот как устроен мир: все, что вы можете сказать, это то, что эти конкретные комбинации, которые выводятся PCA, почти не изменяются в данных. Чтобы справиться с этим, некоторые люди используют самые большие («основные») компоненты непосредственно как независимые переменные в регрессии или последующем анализе, в какой бы форме она ни принималась. Если вы сделаете это, не забудьте сначала удалить зависимую переменную из набора переменных и повторить PCA!

Вот код для воспроизведения этой фигуры:

par(mfrow=c(2,3))
beta <- c(1,1,1,1) # Also can be a matrix with `n.obs` rows: try it!
process(z, beta, sd=0, main="A=B+C+D+E; No error")
process(z, beta, sd=1/10, main="A=B+C+D+E; Small error")
process(z, beta, sd=1/3, threshold=2/3, main="A=B+C+D+E; Large error")

beta <- c(1,1/2,1/2,1)
process(z, beta, sd=0, main="A=B+(C+D)/2+E; No error")
process(z, beta, sd=1/10, main="A=B+(C+D)/2+E; Small error")
process(z, beta, sd=1/3, threshold=2/3, main="A=B+(C+D)/2+E; Large error")

(Мне пришлось поиграться с пороговым значением в случаях с большой ошибкой, чтобы отобразить только один компонент: вот причина для того, чтобы передать это значение в качестве параметра process.)

Пользователь ttnphns любезно направил наше внимание на тесно связанную тему. Один из его ответов (автор JM) предлагает подход, описанный здесь.

Whuber
источник

Вау, вот что я понимаю из твоего ответа ... исправь мои переменные против любой другой переменной. Используйте VIF, чтобы затем идентифицировать связанные переменные .. это работает. Лучше ли делать это с кусками данных одновременно? Также вы удаляете что-нибудь, если вы обнаруживаете колинеарность с использованием предшествующей регрессии? .. Из того, что я понимаю о PCA, как правило, это то, что вы используете самые большие ПК (объясняющие наибольшую дисперсию) на основе собственных значений, так как они объясняют наибольшую дисперсию, они загружаются в различные градусов с использованием исходных переменных. Я не уверен относительно того, что небольшие нагрузки и с чем они коллинеарны

Самуил

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

whuber

"loadings," which are linear combinations of the original variables

A

$\bf A$

A^{- 1}

$\bf A^{-1}$

Кроме того, могу ли я попросить вас оставить свое мнение о возможном использовании операции развертки ( stats.stackexchange.com/a/16391/3277 ) в задаче отслеживания линейно зависимых подмножеств переменных?

ttnphns

X

$X$

X = U W V^{^{'}}

$X=UWV^{'}$

V

$V$ princomp

X V = U W

$XV=UW$

W

$W$

U W

$UW$

0

$0$

X V

$XV$

X

$X$

5

$502\times 480$

JM не является статистиком
источник

3

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

Я создал цикл for (), который вычисляет ранг матрицы по одному столбцу за раз. Таким образом, для первой итерации ранг будет 1. Второй, 2. Это происходит до тех пор, пока ранг не станет МЕНЬШЕ, чем номер столбца, который вы используете.

Очень просто:

for (i in 1:47) {

  print(qr(data.frame[1:i])$rank) 
  print(i) 
  print(colnames(data.frame)[i])
  print("###") 
}

for () разрыв цикла

вычисляет ранг для i-го столбца
печатает номер итерации
печатает имя столбца для справки
делит консоль с "###", чтобы вы могли легко прокрутить

Я уверен, что вы можете добавить оператор if, он мне пока не нужен, потому что я имею дело только с 50-ю столбцами.

Надеюсь это поможет!

Ник П
источник

2

Хотя в этом нет ничего плохого теоретически, это численно нестабильный и неэффективный алгоритм. Особенно с большим количеством столбцов он может не обнаружить почти коллинеарность и ошибочно обнаружить коллинеарность там, где ее нет.

whuber

2

Ранг, r матрицы = количество линейно независимых столбцов (или строк) матрицы. Для п по п матрицы А , ранг (А) = п => все столбцы (или строки) линейно независимы.

Arun
источник

2

Не то чтобы ответ @Whuber действительно нужно было расширить, но я решил дать краткое описание математики.

$\mathbf{X'Xv}=\mathbf{0}$ $\mathbf{v}\neq\mathbf{0}$ $\mathbf{v}$ $\mathbf{X'X}$ $\lambda=0$ $\mathbf{X'X}$ $\mathbf{X}$ $\mathbf{X'X}$ $\lambda=0$ $\mathbf{X'X}$ $\mathbf{v}$ $\lambda$

κ_{j} = \frac{λ_{m a x}}{λ_{j}}

$\mathbf{\kappa_j}=\frac{\lambda_{max}}{\lambda_j}$

X^{'} X = [\begin{array}{rr} 0.001 & 0 & 0 \\ 0 & 0.001 & 0 \\ 0 & 0 & 0.001 \end{array}] .

$\mathbf{X'X}=\left[\begin{array}{rr} 0.001 & 0 & 0 \\ 0 & 0.001 & 0 \\ 0 & 0 & 0.001 \\ \end{array} \right].$

λ_{1} = λ_{2} = λ_{3} = 0.001

$\lambda_1=\lambda_2=\lambda_3=0.001$

κ = \frac{λ_{m a x}}{λ_{m i n}} = 1

$\mathbf{\kappa}=\frac{\lambda_{max}}{\lambda_{min}}=1$

Цитирование

Монтгомери, Д. (2012). Введение в анализ линейной регрессии, 5-е издание. John Wiley & Sons Inc.

tjnel
источник

1

X

$X$

Q R

$QR$

X

$X$

n

$n$

k

$k$

n >> k

$n >> k$

X

$X$

R

$R$

Тестирование на линейную зависимость среди столбцов матрицы

Ответы: