python - Cómo encontrar filas/columnas degeneradas en una matriz de covarianza

matrix numpy

Question

Estoy usando numpy.cov para crear una matriz de covarianza a partir de un conjunto de datos de más de 400 series temporales. El uso de linalg.det me da un valor de cero, por lo que la matriz es singular. Puedo usar linalg.svd para ver que el rango es dos menos que el número de columnas, así que en algún lugar de la matriz de covarianza tengo algunas combinaciones lineales para hacer que la matriz se degenere. He usado corrcoef en las series temporales subyacentes, pero no hay correlación> 0,78, por lo que no es obvio. ¿Alguien puede sugerir un método para determinar la ubicación de las columnas degeneradas? Gracias.

Answer 1

Si toma la descomposición QR de una matriz A , las columnas de R con un valor distinto de cero a lo largo de la diagonal corresponden a columnas linealmente independientes de A

import numpy as np
linalg = np.linalg

def independent_columns(A, tol = 1e-05):
    """
    Return an array composed of independent columns of A.

    Note the answer may not be unique; this function returns one of many
    possible answers.

    http://stackoverflow.com/q/13312498/190597 (user1812712)
    http://math.stackexchange.com/a/199132/1140 (Gerry Myerson)
    http://mail.scipy.org/pipermail/numpy-discussion/2008-November/038705.html
        (Anne Archibald)

    >>> A = np.array([(2,4,1,3),(-1,-2,1,0),(0,0,2,2),(3,6,2,5)])
    >>> independent_columns(A)
    np.array([[1, 4],
              [2, 5],
              [3, 6]])
    """
    Q, R = linalg.qr(A)
    independent = np.where(np.abs(R.diagonal()) > tol)[0]
    return A[:, independent]

def matrixrank(A,tol=1e-8):
    """
    http://mail.scipy.org/pipermail/numpy-discussion/2008-February/031218.html
    """
    s = linalg.svd(A,compute_uv=0)
    return sum( np.where( s>tol, 1, 0 ) )


matrices = [
    np.array([(2,4,1,3),(-1,-2,1,0),(0,0,2,2),(3,6,2,5)]),
    np.array([(1,2,3),(2,4,6),(4,5,6)]).T,
    np.array([(1,2,3,1),(2,4,6,2),(4,5,6,3)]).T,
    np.array([(1,2,3,1),(2,4,6,3),(4,5,6,3)]).T,
    np.array([(1,2,3),(2,4,6),(4,5,6),(7,8,9)]).T
    ]

for A in matrices:
    B = independent_columns(A)
    assert matrixrank(A) == matrixrank(B) == B.shape[-1]

assert matrixrank(A) == matrixrank(B) comprueba que la función independent_columns devuelve una matriz del mismo rango que A

assert matrixrank(B) == B.shape[-1] verifica que el rango de B igual al número de columnas de B