2012-10-26 9 views
8

Tôi có ma trận mxn lớn và tôi đã xác định các cột phụ thuộc tuyến tính. Tuy nhiên, tôi muốn biết nếu có một cách trong R để viết các cột phụ thuộc tuyến tính về những người độc lập tuyến tính. Vì nó là một ma trận lớn nên không thể thực hiện dựa trên việc kiểm tra.Làm cách nào để viết cột phụ thuộc tuyến tính trong ma trận theo các cột độc lập tuyến tính?

Đây là ví dụ về loại ma trận mà tôi có.

> 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 

Đây rõ ràng là x3 = x1-x2, x5 = x1-x4. Tôi muốn biết nếu có một cách tự động để có được điều đó cho một ma trận lớn hơn.

Cảm ơn!

+1

Chức năng này có thể giúp bạn: http://www.inside-r.org/packages/cran/heplots/docs/gsorth –

+1

@BenBolker Và kể từ khi liên kết đó là bây giờ đã chết chỉ cần tìm hàm 'gsorth' ở đây: https://cran.r-project.org/web/packages/heplots/heplots.pdf – Dason

Trả lời

9

Tôi chắc chắn có một cách tốt hơn nhưng tôi cảm thấy thích chơi đùa với điều này. Tôi về cơ bản làm một kiểm tra ở đầu để xem nếu ma trận đầu vào là cột xếp hạng đầy đủ để tránh tính toán không cần thiết trong trường hợp nó được xếp hạng đầy đủ. Sau đó tôi bắt đầu với hai cột đầu tiên và kiểm tra xem submatrix đó có xếp hạng cột đầy đủ hay không, nếu đó là sau đó tôi kiểm tra các cột đầu tiên và vv. Một khi chúng ta tìm thấy một số submatrix không có thứ hạng cột đầy đủ, tôi sẽ lấy lại cột cuối cùng trong submatrix đó trong phần trước, nó cho chúng ta biết cách xây dựng tổ hợp tuyến tính của các cột đầu tiên để lấy cột cuối cùng.

Chức năng của tôi hiện không sạch sẽ và có thể thực hiện một số kiểm tra bổ sung nhưng ít nhất đó là một sự khởi đầu.

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) 

mang đến cho

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

Cảm ơn !!! Điều đó hoạt động tuyệt vời! –

+0

@ Dưa tôi đã đi qua câu trả lời này trong khi nghiên cứu cho câu hỏi của tôi (rất giống, chỉ với các yêu cầu bổ sung) @ http://stackoverflow.com/questions/43596486/how-to-identify-columns-that-are-sums- of-other-columns-in-a-dataset. Sẽ rất tuyệt vời khi bạn có trọng số trong cuộc thảo luận đó. – Aurimas