Test de dépendance linéaire entre les colonnes d'une matrice

26

J'ai une matrice de corrélation des retours de titres dont le déterminant est zéro. (Cela est un peu surprenant car la matrice de corrélation d'échantillon et la matrice de covariance correspondante devraient théoriquement être définies positives.)

Mon hypothèse est qu'au moins un titre dépend linéairement d'autres titres. Y a-t-il une fonction dans R qui teste séquentiellement chaque colonne une matrice pour la dépendance linéaire?

Par exemple, une approche consisterait à constituer une matrice de corrélation un titre à la fois et à calculer le déterminant à chaque étape. Lorsque le déterminant = 0, arrêtez-vous car vous avez identifié le titre qui est une combinaison linéaire d'autres titres.

Toute autre technique permettant d'identifier la dépendance linéaire dans une telle matrice est appréciée.

Ram Ahluwalia
la source
Votre matrice est semi-définie positive, mais elle n'est pas définie positive, car elle est singulière.
ttnphns
Quelles sont les dimensions (no. Variables; no. Samples)?
Karl
Nombre de colonnes = 480. Nombre de lignes pour chaque série chronologique = 502. En général, vous constatez que plus la série chronologique est grande, plus la matrice de covariance de l'échantillon a tendance à être définie positive. Cependant, il existe de nombreux cas où vous souhaitez utiliser une valeur de T beaucoup plus petite (ou un poids exponentiel) pour refléter les conditions récentes du marché.
Ram Ahluwalia du
3
La question est mal posée. Si votre matrice de données est 480 par 502, dire que la matrice a un rang (l'espace de colonne de la matrice a la dimension ) équivaut mathématiquement à dire qu'une colonne est une combinaison linéaire des autres, mais vous pouvez 'choisissez pas une colonne et dites que c'est la colonne qui dépend linéairement. Il n'y a donc pas de procédure pour ce faire, et la procédure suggérée choisira une sécurité assez arbitraire en fonction de l'ordre dans lequel ils sont inclus. q < 480q<480q<480
NRH
La matrice de covariance est symétrique. Il est généré par transposition (A) * A. La matrice A a des dimensions 480x502. Cependant, la matrice de covariance est 480x480
Ram Ahluwalia

Réponses:

6

Vous semblez poser une question vraiment provocante: comment détecter, étant donné une matrice de corrélation (ou covariance, ou somme des carrés et produits croisés) singulière, quelle colonne dépend linéairement de laquelle. Je suppose que l' opération de balayage pourrait aider. Voici ma sonde en SPSS (pas R) pour illustrer.

Générons quelques données:

        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

Créons une dépendance linéaire entre V2, V4 et V5:

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

Nous avons donc modifié notre colonne 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.

Les impressions de M en 5 itérations:

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

Remarquez que finalement la colonne 5 est pleine de zéros. Cela signifie (si je comprends bien) que V5 est lié linéairement avec certaines des colonnes précédentes . Quelles colonnes? Regardez l'itération où la colonne 5 n'est pas remplie de zéros en dernier - itération 4. Nous voyons là que V5 est lié avec V2 et V4 avec des coefficients -.3333 et .8333: V5 = -.3333 * V2 + .8333 * V4, ce qui correspond à ce que nous avons fait avec les données: V4 = .4 * V2 + 1.2 * V5.

C'est ainsi que nous avons su quelle colonne est liée linéairement avec quelle autre. Je n'ai pas vérifié l'utilité de l'approche ci-dessus dans un cas plus général avec de nombreux groupes d'interdépendances dans les données. Dans l'exemple ci-dessus, cela semble cependant utile.

ttnphns
la source
N'est-ce pas la forme échelonnée de la ligne réduite? Si oui, n'y a-t-il pas de packages / fonctions disponibles dans R?
Arun
@Arun, je ne suis pas un utilisateur R, je ne peux donc pas le savoir.
ttnphns
25

Voici une approche simple: calculer le rang de la matrice qui résulte de la suppression de chacune des colonnes. Les colonnes qui, lorsqu'elles sont supprimées, donnent le rang le plus élevé sont celles qui dépendent linéairement (puisque la suppression de celles-ci ne diminue pas le rang, tandis que la suppression d'une colonne linéairement indépendante le fait).

Dans R:

rankifremoved <- sapply(1:ncol(your.matrix), function (x) qr(your.matrix[,-x])$rank)
which(rankifremoved == max(rankifremoved))
James
la source
1
Réponse extrêmement utile pour déterminer la colonne incriminée dans une matrice de régression où j'ai reçu l'erreur system is exactly singular: U[5,5] = 0 , ce que je sais maintenant signifie que la colonne 5 était le problème (semble évident avec le recul car c'est une colonne de zéros!)
Matt Weller
Dans le commentaire de James, il a posté le script: rankifremoved <- sapply (1: ncol (your.matrix), fonction (x) qr (your.matrix [, - x]) $ rank) qui (rankifremoved == max ( rankifremoved)) J'ai fait un test sur une matrice, je voudrais savoir sur la sortie de R. Les colonnes de la sortie sont linéairement dépendantes? Remercier!
@ EltonAraújo: La sortie sera un vecteur donnant les indices des colonnes linéairement dépendantes: donc (2,4,5) pour l'exemple dans la réponse de ttnphns. Mais je me demande comment les problèmes de précision numérique vont affecter cette méthode.
Scortchi - Réintégrer Monica
rankifremoved contient toutes les colonnes qui dépendent linéairement entre elles ou entre elles. Dans certaines applications, nous pouvons vouloir conserver une colonne ou quelques colonnes et ne pas tout supprimer
MasterJedi
Cela ne devrait-il pas renvoyer un ensemble vide your.matrix = matrix(1:4, 2)?
Holger Brandl
15

La question porte sur «l'identification des relations [linéaires] sous-jacentes» entre les variables.

Le moyen rapide et facile de détecter des relations consiste à régresser toute autre variable (utiliser une constante, même) par rapport à ces variables à l'aide de votre logiciel préféré: toute bonne procédure de régression détectera et diagnostiquera la colinéarité. (Vous n'aurez même pas la peine de regarder les résultats de la régression: nous comptons simplement sur un effet secondaire utile de la configuration et de l'analyse de la matrice de régression.)

0

(Il y a un art et beaucoup de littérature associé à l'identification de ce qu'est un "petit" chargement. Pour modéliser une variable dépendante, je suggère de l'inclure dans les variables indépendantes de l'ACP afin d'identifier les composants - indépendamment de leur taille - dans laquelle la variable dépendante joue un rôle important. De ce point de vue, "petit" signifie beaucoup plus petit que n'importe quel composant de ce type.)


Regardons quelques exemples. (Celles-ci sont utilisées Rpour les calculs et le traçage.) Commencez par une fonction pour effectuer l'ACP, recherchez les petits composants, tracez-les et retournez les relations linéaires entre eux.

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,,EUNE

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

B,,EUNE=B+C++EUNE=B+(C+)/2+Esweep

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,,EUNE

Résultats

La sortie associée au panneau supérieur gauche était

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

00UNE-B-C--E

La sortie du panneau central supérieur était

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

(UNE,B,C,,E)

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

UNE=B+C++E

1,1/2,1/2,1

En pratique, il n’est souvent pas vrai qu’une variable soit désignée comme une combinaison évidente des autres: tous les coefficients peuvent être de tailles comparables et de signes variables. De plus, lorsqu'il existe plusieurs dimensions de relations, il n'y a pas de moyen unique de les spécifier: une analyse plus approfondie (comme la réduction des lignes) est nécessaire pour identifier une base utile pour ces relations. C'est ainsi que le monde fonctionne: tout ce que vous pouvez dire, c'est que ces combinaisons particulières qui sont produites par PCA ne correspondent à presque aucune variation dans les données. Pour y faire face, certaines personnes utilisent directement les plus grandes composantes («principales») comme variables indépendantes dans la régression ou l'analyse subséquente, quelle que soit la forme qu'elle pourrait prendre. Si vous faites cela, n'oubliez pas d'abord de supprimer la variable dépendante de l'ensemble de variables et de refaire le PCA!


Voici le code pour reproduire cette figure:

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")

(J'ai dû jouer avec le seuil dans les cas de grandes erreurs afin d'afficher un seul composant: c'est la raison pour laquelle je fournis cette valeur comme paramètre process.)


L'utilisateur ttnphns a aimablement dirigé notre attention sur un sujet étroitement lié. Une de ses réponses (par JM) suggère l'approche décrite ici.

whuber
la source
Wow, voici ce que je comprends de votre réponse .. régresser mes variables contre toute autre variable. Utilisez le VIF pour ensuite identifier les variables liées ... cela fonctionne. Est-il préférable de le faire avec des morceaux de données à la fois? Supprimez-vous également quoi que ce soit si vous détectez une colinéarité en utilisant la régression avant? degrés en utilisant les variables d'origine. Je ne sais pas avec quels petits chargements et avec quoi ils sont colinéaires
Samuel
Cette réponse explique comment interpréter les petits composants: ils présentent les colinéarités. Oui, vous pouvez utiliser des sous-groupes de variables si vous le souhaitez. La méthode de régression consiste simplement à détecter la présence de colinéarité, non à identifier les relations colinéaires: c'est ce que fait l'ACP.
whuber
"loadings," which are linear combinations of the original variablesUNEUNE-1
De plus, puis-je vous demander de laisser votre avis sur l'utilisation possible de l'opération de balayage ( stats.stackexchange.com/a/16391/3277 ) dans la tâche de traquer des sous-ensembles de variables linéairement dépendants?
ttnphns
XX=UWVVprincompXV=UWWUW0XVX
5

502×480

JM n'est pas un statisticien
la source
3

J'ai rencontré ce problème il y a environ deux semaines et j'ai décidé que je devais le réexaminer car, lorsqu'il s'agit d'ensembles de données massifs, il est impossible de faire ces choses manuellement.

J'ai créé une boucle for () qui calcule le rang de la matrice une colonne à la fois. Ainsi, pour la première itération, le classement sera 1. Le second, 2. Cela se produit jusqu'à ce que le classement devienne MOINS que le numéro de colonne que vous utilisez.

Très simple:

for (i in 1:47) {

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

pour () rupture de boucle

  1. calcule le rang de la ième colonne
  2. imprime le numéro d'itération
  3. imprime le nom de la colonne pour référence
  4. divise la console avec "###" afin que vous puissiez facilement faire défiler

Je suis sûr que vous pouvez ajouter une instruction if, je n'en ai pas encore besoin car je ne traite que de 50 colonnes.

J'espère que cela t'aides!

Nick P
la source
2
Bien qu'il n'y ait rien de mal à cela théoriquement, c'est un algorithme numériquement instable et inefficace. Surtout avec un grand nombre de colonnes, il peut ne pas détecter la quasi-colinéarité et détecter faussement la colinéarité là où il n'en existe pas.
whuber
2

Rang, r d'une matrice = nombre de colonnes (ou lignes) linéairement indépendantes d'une matrice. Pour une matrice n par n A , le rang (A) = n => toutes les colonnes (ou lignes) sont linéairement indépendantes.

Arun
la source
2

Non pas que la réponse donnée par @Whuber doive vraiment être développée, mais j'ai pensé fournir une brève description des mathématiques.

XXv=0v0vXXλ=0XXXXXλ=0XXvλ

κj=λmuneXλj

XX=[0,0010000,0010000,001].
λ1=λ2=λ3=0,001
κ=λmuneXλmjen=1

Citations

Montgomery, D. (2012). Introduction à l'analyse de régression linéaire, 5e édition. John Wiley & Sons Inc.

tjnel
la source
1
X
QRXnkn>>kXR