Got Questions? Get Answers.
Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
well conditioned basis of rectangular matrix A

Subject: well conditioned basis of rectangular matrix A

From: HenryW

Date: 25 Mar, 2009 12:38:57

Message: 1 of 7

I have an m by n matrix A, m << n. It is full row rank and I
would like to find a set of n columns of A that form a basis
of the range of A so that the columns (in matrix B say) are well conditioned, i.e. B is well conditioned. I know that there is code out there, e.g. MA48, that does this. But is there some similar code in matlab?

Subject: well conditioned basis of rectangular matrix A

From: Bruno Luong

Date: 25 Mar, 2009 13:21:02

Message: 2 of 7

HenryW <hwolkowicz@uwaterloo.ca> wrote in message <29402152.1237984780487.JavaMail.jakarta@nitrogen.mathforum.org>...
> I have an m by n matrix A, m << n. It is full row rank and I
> would like to find a set of n columns of A that form a basis
> of the range of A so that the columns (in matrix B say) are well conditioned, i.e. B is well conditioned. I know that there is code out there, e.g. MA48, that does this. But is there some similar code in matlab?

Interesting problem. What are typical values of m and n?

Do you need the optimal solution (maximizing the condition number) or just sub-optimal?

Could you explain what it will be used for?

Bruno

Subject: well conditioned basis of rectangular matrix A

From: John D'Errico

Date: 25 Mar, 2009 13:52:01

Message: 3 of 7

HenryW <hwolkowicz@uwaterloo.ca> wrote in message <29402152.1237984780487.JavaMail.jakarta@nitrogen.mathforum.org>...
> I have an m by n matrix A, m << n. It is full row rank and I
> would like to find a set of n columns of A that form a basis
> of the range of A so that the columns (in matrix B say) are well conditioned, i.e. B is well conditioned. I know that there is code out there, e.g. MA48, that does this. But is there some similar code in matlab?

orth does it simply, if you are willing to use an
othogonal basis, so not the original columns of A.

If you just want a set of the best columns, then
use a QR with column pivoting, throwing away
Q and R. All you need are the column pivots to
identify a set of columns that are best for this
purpose. A problem will arise when you need to
decide when to stop, but looking at the diagonal
elements of R will help you to make that decision
too. Stop when they get as small as you wish.

A = rand(10,5);
A = [A,A+1];
A = A(:,randperm(10));

As built, A will have rank 6. Can we find a
subset of the columns of A that are maximally
independent?

[Q,R,E] = qr(A,0);

E =
     8 4 3 2 6 10 9 1 7 5

diag(R)
ans =
       -4.898
      -1.4934
       -1.228
      0.87374
      0.59921
      0.24203
  -8.8239e-16
  -2.8686e-16
   -2.322e-16
   1.6598e-16

As you can see from the diagonal elements of R,
the rank of A is 6. This agrees with rank.

rank(A)
ans =
     6

Just choose the columns given by the first 6
elements of E for a quite reasonable choice.

svd(A(:,E(1:6)))
ans =
       10.395
         1.37
       1.0433
      0.88079
        0.523
      0.16629

It will be more difficult if A does not have so
clear a numerical rank. But then that just makes
life interesting.

HTH,
John

Subject: well conditioned basis of rectangular matrix A

From: Bruno Luong

Date: 25 Mar, 2009 14:04:01

Message: 4 of 7

But John,

Does qr gives what OP really wants?

For example:
 
A=[1 0 1e-6 0;
     0 1e6 0 1e-6]

The selection based on QR decomposition will pick column [1 2]
and the condition number is 1e6.

If we select column [3 4], the condition number is 1.

Bruno

Subject: well conditioned basis of rectangular matrix A

From: John D'Errico

Date: 25 Mar, 2009 14:27:01

Message: 5 of 7

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <gqddkh$8gh$1@fred.mathworks.com>...
> But John,
>
> Does qr gives what OP really wants?
>
> For example:
>
> A=[1 0 1e-6 0;
> 0 1e6 0 1e-6]
>
> The selection based on QR decomposition will pick column [1 2]
> and the condition number is 1e6.
>
> If we select column [3 4], the condition number is 1.
>
> Bruno

But knowing what a poster REALLY wants is so
difficult. ;-)

Are those the best columns to choose? In some
respects they might be. In other ways, not so.
It very much depends upon his goals and what
he knows about his matrix and what he will do
with the identified columns. Perhaps an initial
column scaling so all have a unit norm might
be appropriate.

John

Subject: well conditioned basis of rectangular matrix A

From: HenryW

Date: 25 Mar, 2009 14:25:42

Message: 6 of 7

m is typically of the order 10^2 and n is of the order 10^5
This comes from semidefinite programming, SDP, constraints
         A X = b
X here is a symmetric matrix and A the linear transformation into R^m.
The factorization is used to get a representation
       X = Xhat + V y,
so that we can replace the X >=0 constraint with
    -Vy <= Xhat
This allows us to work with the 'dual formulation' of the primal problem. (We do not want to work with the dual of the problem since in SDP there can be duality gaps.)

Subject: well conditioned basis of rectangular matrix A

From: HenryW

Date: 25 Mar, 2009 14:49:49

Message: 7 of 7

I received the following reply from a private email - so I cannot take the credit.

for finding a good partition A = (B E)
where B is a basis as well-conditioned as A allows,
an efficient and reliable approach is to use a sparse LU
factorization of A(transpose) in which L is well-conditioned.

In Matlab, this means using the *original* sparse LU as follows:


   [m,n] = size(A);
   AT = A';
   q = colamd(AT);
   [L,U,p] = lu(AT(:,q),'vector');
   B = A(:,p(1:m));
   E = A(:,p(m+1:n));


Beware that [L,U,P,Q] = lu(AT,...) does not guarantee a
well-conditioned L and cannot be used in this context.
-----------------------------------------
In terms of what I need - I have an optimization problem with a constraints
   Ax=b x>=0
and would like to replace it with
      inv(B) N y <= inv(B) b, x_N >=0
so I want B well conditioned so that
inv(B) N is easily found and accurate.
I do not hope for sparsity here though.
    thanks!!

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us