Come scrivere una colonna dipendente linearmente in una matrice in termini di colonne linearmente indipendenti?

Question

Ho una grande matrice mxn e ho identificato le colonne dipendenti linearmente. Tuttavia, voglio sapere se c'è un modo in R per scrivere le colonne dipendenti linearmente in termini di quelle linearmente indipendenti. Dal momento che è una matrice di grandi dimensioni, non è possibile fare in base all'ispezione.

Ecco un esempio di giocattolo del tipo di matrice che ho.

> mat <- matrix(c(1,1,0,1,0,1,1,0,0,1,1,0,1,1,0,1,0,1,0,1), byrow=TRUE, ncol=5, nrow=4) 
> mat 
     [,1] [,2] [,3] [,4] [,5] 
[1,] 1 1 0 1 0 
[2,] 1 1 0 0 1 
[3,] 1 0 1 1 0 
[4,] 1 0 1 0 1 

Qui è evidente che x3 = x1-x2, x5 = x1-x4. Voglio sapere se c'è un modo automatico per ottenerlo per una matrice più grande.

Grazie!

Answer 1

Sono sicuro che c'è un modo migliore ma mi sembrava di giocare con questo. Fondamentalmente faccio un controllo all'inizio per vedere se la matrice di input è piena di colonne per evitare calcoli inutili nel caso in cui sia completa. Dopo di ciò, comincio con le prime due colonne e controllo se quella submatrix è di rango completo di colonne, se è allora controllo le prime colonne e così via. Una volta che abbiamo trovato qualche submatrix che non ha il rank di colonne complete, regredisco l'ultima colonna di quella submatrix sulla precedente che ci dice come costruire combinazioni lineari delle prime colonne per ottenere l'ultima colonna.

La mia funzione non è molto pulita in questo momento e potrebbe fare qualche controllo aggiuntivo ma almeno è un inizio.

mat <- matrix(c(1,1,0,1,0,1,1,0,0,1,1,0,1,1,0,1,0,1,0,1), byrow=TRUE, ncol=5, nrow=4) 


linfinder <- function(mat){ 
    # If the matrix is full rank then we're done 
    if(qr(mat)$rank == ncol(mat)){ 
     print("Matrix is of full rank") 
     return(invisible(seq(ncol(mat)))) 
    } 
    m <- ncol(mat) 
    # cols keeps track of which columns are linearly independent 
    cols <- 1 
    for(i in seq(2, m)){ 
     ids <- c(cols, i) 
     mymat <- mat[, ids] 
     if(qr(mymat)$rank != length(ids)){ 
      # Regression the column of interest on the previous 
      # columns to figure out the relationship 
      o <- lm(mat[,i] ~ mat[,cols] + 0) 
      # Construct the output message 
      start <- paste0("Column_", i, " = ") 
      # Which coefs are nonzero 
      nz <- !(abs(coef(o)) <= .Machine$double.eps^0.5) 
      tmp <- paste("Column", cols[nz], sep = "_") 
      vals <- paste(coef(o)[nz], tmp, sep = "*", collapse = " + ") 
      message <- paste0(start, vals) 
      print(message) 
     }else{ 
      # If the matrix subset was of full rank 
      # then the newest column in linearly independent 
      # so add it to the cols list 
      cols <- ids 
     } 
    } 
    return(invisible(cols)) 
} 

linfinder(mat) 

che dà

> linfinder(mat) 
[1] "Column_3 = 1*Column_1 + -1*Column_2" 
[1] "Column_5 = 1*Column_1 + -1*Column_4"