0

I am trying to find independent columns to solve the system of linear equations. Here my simplified example:

> mat = matrix(c(1,0,0,0,-1,1,0,0,0,-1,1,0,0,0,-1,0,-1,0,0,1,0,0,1,-1), nrow=4, ncol=6, dimnames=list(c("A", "B", "C", "D"), paste("v", 1:6, sep="")))
> mat
  v1 v2 v3 v4 v5 v6
A  1 -1  0  0 -1  0
B  0  1 -1  0  0  0
C  0  0  1 -1  0  1
D  0  0  0  0  1 -1

The matrix is full rank:

qr(mat)$rank

gives me 4, and since there are 6 columns, there should be 6-4=2 independent columns from which I can calculate the others. I know that columns v4 and v6 are independent... My first question is, how can I find these columns (maybe with qr(mat)$pivot)?

By rearranging the linear equations on paper, I see that
[v1, v2, v3, v4, v5, v6] = [v4, v4-v6, v4-v6, v4, v4, v6, v6]

and thus I can find from arbitrary values for v4 and v6 a vector that lies in the null space by multiplying v4 and v6 with the vectors below:

v4 * [1,1,1,1,0,0] + v6 * [0,-1,-1,0,1,1]

My second question is: How do I find these vectors, meaning how do I solve the matrix for v4 and v6? For example

qr.solve(mat, cbind(c(0,0,0,0), c(0,0,0,0)))

gives me two vectors of length 6 with only zeros.

Any help is appreciated, many thanks in advance!

-H-

user000001
  • 32,226
  • 12
  • 81
  • 108
user1981275
  • 13,002
  • 8
  • 72
  • 101
  • No. Since rank is 4 there are 4 independent columns. Furthermore, it's not as though 2 specific ones are dependent, only that if you pick 3 of them then only one more can be picked that will be also independent. Unless there are a pair that are simple multiples, then you might be able to use any one of them as a basis vector. – IRTFM Feb 18 '13 at 19:00
  • Are then v1, v2, v3, v5 the independent columns? But I can solve the system of equations with only v4 and v6 as independent variables, can't I? – user1981275 Feb 18 '13 at 19:04
  • They are, but you could also have picked v1, v2, v4, v6. – IRTFM Feb 18 '13 at 19:07
  • OK, I just did id on paper and indeed I can take v1, v2, v4 and v6 and form basis vectors with v3 and v5. I picked v4 and v6 because these columns did not contain pivots. – user1981275 Feb 18 '13 at 19:18

2 Answers2

3

Use the pivot information to find a set of independent columns:

q <- qr(mat)

mmat <- mat[,q$pivot[seq(q$rank)]]

mmat
##   v1 v2 v3 v5
## A  1 -1  0 -1
## B  0  1 -1  0
## C  0  0  1  0
## D  0  0  0  1

qr(mmat)$rank
## [1] 4

Why does this work? The meaning of pivot is given in QR.Auxiliaries {base} brought up with ?qr.Q. In particular:

qr.R returns R. This may be pivoted, e.g., if a <- qr(x) then x[, a$pivot] = QR.
The number of rows of R is either nrow(X) or ncol(X) (and may depend on whether
complete is TRUE or FALSE).

Pivoting is done to order the eigenvalues in decreasing absolute value, for numerical stability. This also means that any 0 eigenvalues are at the end, beyond q$rank in q$pivot (and nonexistent in the current example, where Q is a 4x4 orthogonal matrix).

The final lines in the QR.Auxiliaries {base} show this relationship:

pivI <- sort.list(a$pivot) # the inverse permutation
stopifnot(
 all.equal(x[, a$pivot], qr.Q(a) %*% qr.R(a)),          # TRUE
 all.equal(x           , qr.Q(a) %*% qr.R(a)[, pivI]))  # TRUE too!
Matthew Lundberg
  • 42,009
  • 6
  • 90
  • 112
  • This seems to work fine, although I don't yet fully understand what q$pivot actually gives me, ?qr says it is "information on the pivoting strategy used during the decomposition" – user1981275 Feb 18 '13 at 19:20
  • Now to find the vectors [1,1,1,1,0,0] for v4 and [0,-1,-1,0,1,1] for v6, do I have to solve mmat? – user1981275 Feb 18 '13 at 19:22
  • Where exactly can I verify that pivoting orders the eigenvalues in decreasing absolute value? – landau Apr 13 '22 at 12:27
0

If you start with v4 and v6 then you need 2 more with non-zero values inrows 1 and 2 so that you need to pick v1 and either v2 or v3. These are all possible basis choices that will have maximal rank.

> qr(mat[, c(1,2,4,6)])$rank
[1] 4
> qr(mat[, c(1,2,3,5)])$rank
[1] 4
> qr(mat[, c(1,3,4,6)])$rank
[1] 4

It is simply not the case that "independent columns" are uniquely determined. There may be sets of columns that are necessarily dependent, e.g., ones which are scalar multiples of each other, but that is not the case here.

On the other hand this will be rank deficient:

> qr(mat[, c(1,2,3,4)])$rank
[1] 3
IRTFM
  • 258,963
  • 21
  • 364
  • 487